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Abstract 



The study of matter at extreme densities has been a major focus in theoretical physics in 
the last half-century. The wide spectrum of information that the field produces provides an 
invaluable contribution to our knowledge of the world in which we live. Most fascinatingly, 
the insight into the world around us is provided from knowledge of the intangible, at both 
the smallest and largest scales in existence. 

Through the study of nuclear physics we are able to investigate the fundamental con- 
struction of individual particles forming nuclei, and with further physics we can extrapolate 
to neutron stars. The models and concepts put forward by the study of nuclear matter help 
to solve the mystery of the most powerful interaction in the universe; the strong force. 

In this study we have investigated a particular state-of-the-art model which is currently 
used to refine our knowledge of the workings of the strong interaction and the way that it 
is manifested in both neutron stars and heavy nuclei, although we have placed emphasis on 
the former for reasons of personal interest. The main body of this work has surrounded an 
effective field theory known as Quantum Hadrodynamics (QHD) and its variations, as well 
as an extension to this known as the Quark-Meson Coupling (QMC) model, and variations 
thereof. We further extend these frameworks to include the possibility of a phase transition 
from hadronic matter to deconfined quark matter to produce hybrid stars, using various 
models. 

We have investigated these pre-existing models to deeply understand how they are jus- 
tified, and given this information, we have expanded them to incorporate a modern under- 
standing of how the strong interaction is manifest. 
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Introduction 



The world we live in is a strange and wonderful place. The vast majority of our interactions 
with objects around us can be described and predicted by relatively simple equations 
and relations, which have been fully understood for centuries now. The concept of this tactile 
world is however, a result of averaging a vast ensemble of smaller effects on an unimaginably 
small scale, all conspiring together to produce what we observe at the macroscopic level. In 
order to investigate this realm in which particles can fluctuate in and out of existence in the 
blink of an eye, we must turn to more sophisticated descriptions of the players on this stage. 

This research has the goal of furthering our understanding of the interior of what are 
loosely called 'neutron stars', though as shall be shown herein, the contents are not necessarily 
only neutrons. With this in mind, further definitions include ^hyperon stars' , ^quark stars\ 
and ^hybrid stars\ to describe compact stellar objects containing hyperons^, quarks, and a 
mixture of each, respectively. In order to do this, we require physics beyond that which 
describes the interactions of our daily lives; we require physics that describes the individual 
interactions between particles, and physics that describes the interactions between enormous 
quantities of particles. 

Only by uniting the physics describing the realm of the large and that of the small can one 
contemplate so many orders of magnitude in scale; from individual particles with a diameter 
of less than 10~^^ m, up to neutron stars with a diameter of tens of kilometers. Yet the 
physics at each end of this massive scale are unified in this field of nuclear matter in which 
interactions of the smallest theorized entities conspire in such a way that densities equivalent 
to the mass of humanity compressed to the size of a mere sugar cube become commonplace, 
energetically favourable, and stable. 

The sophistication of the physics used to describe the world of the tiny and that of 
the enormous has seen much development over time. The current knowledge of particles 
has reached a point where we are able to make incredibly precise predictions about the 
properties of single particles and have them confirmed with equally astonishing accuracy 
from experiments. The physics describing neutron stars has progressed from relatively simple 
(yet sufficiently consistent with experiment) descriptions of neutron (and nucleon) matter to 
many more sophisticated descriptions involving various species of baryons, mesons, leptons 
and even quarks. 

The outcome of work such as this is hopefully a better understanding of matter at both 
the microscopic and macroscopic scales, as well as the theory and formalism that unites 
these two extremes. The primary methods which we have used to construct models in this 
thesis are Quantum Hadrodynamics (QHD) — which shall be described in Sec. 3.1 — and the 
Quark-Meson Coupling (QMC) model, which shall be described in Sec. 3.2. 

In this thesis, we will outline the research undertaken in which we produce a model for 
neutron star structure which complies with current theories for dense matter at and above 
nuclear density and is consistent with current data for both finite nuclei and observed neutron 
stars. Although only experimental evidence can successfully validate any theory, we hope 

^Baryons for which one or more of the three valence quarks is a strange quark. For example, the 
hyperon contains two up quarks and a strange quark. For a table of particle properties, including quark 
content, refer to Table B.l 
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1.1. The Four Forces 



to convey a framework and model that possesses a minimum content of inconsistencies or 
unjustified assumptions such that any predictions that are later shown to be fallacious can 
only be attributed to incorrect initial conditions. 

As a final defense of any inconsistencies that may arise between this research and exper- 
iment, we refer the reader to one of the author's favourite quotes: 

There is a theory which states that if ever anyone discovers exactly what 

the Universe is for and why it is here, it will instantly disappear and be replaced 
by somethins even more bizarre and inexplicable. 

There is another theory which states that this has already happened. " 

- Douglas Adams, The Hitchhiker's Guide To The Galaxy. 

Our calculations will begin at the particle interaction level, which will be described in 
more detail in Section 2, from which we are able to reproduce the bulk properties of matter 
at high densities, and which we shall discuss in Sections 3.1-3.5. The methods for producing 
our simulations of the interactions and bulk properties will be detailed in Section 4, along with 
a discussion on how this is applied to the study of compact stellar objects. The results of the 
simulations and calculations will be discussed in detail in Section 5, followed by discussions on 
the interpretation of these results in Section 6. For the convenience of the reader, and for the 
sake of completeness, derivations for the majority of the equations used herein are provided 
in Appendix A, and useful information regarding particles is provided in Appendix B. For 
now however, we will provide a brief introduction to this field of study. 

1.1 The Four Forces 

Theoretical particle physics has seen much success and found many useful applications; from 
calculating the individual properties of particles to precisions that rival even the best ex- 
perimental setups, to determining the properties of ensembles of particles of greater and 
greater scale, and eventually to the properties of macroscopic objects as described by their 
constituents. 

In order to do this, we need to understand each of the four fundamental forces in Nature. 
The weakest of these forces — gravity — attracts any two masses, and will become most im- 
portant in the following section. Slightly stronger is electromagnetism; the force responsible 
for electric charge and magnetism. This force provides an attraction between opposite electric 
charges (and of course, repulsion between like charges), and thus helps to bind electrons to 
nuclei. The mathematical description of this effect is Quantum Electrodynamics (QED). 

The 'weak nuclear force' — often abbreviated to 'the weak force' — is responsible for the 
decay of particles and thus radioactivity in general. At high enough energies, this force is 
unified with electromagnetism into the 'electro- weak force'. The strongest of all the forces 
is the 'strong nuclear force', abbreviated to 'the strong force'. This force is responsible for 
attraction between certain individual particles over a very short scale, and is responsible for 
the binding of protons within nuclei which would otherwise be thrown apart by the repulsive 
electromagnetic force between the positively charged protons. Each of these forces plays a 
part in the work contained herein, but the focus of our study will be the strong force. 
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A remarkably successful description of the strongest force at the microscopic level is 
'Quantum Chromodynamics' (QCD) which is widely believed to be the true description of 
strong interactions, relying on quark and gluon degrees of freedom. The major challenge 
of this theory is that at low energies it is non-perturbative^ , in that the coupling constant 
— which one would normally perform a series-expansion in powers of — is large, and thus is 
not suitable for such an expansion. Regularisation techniques have been produced to create 
perturbative descriptions of QCD, but perturbative techniques fail to describe both dynamical 
chiral symmetry breaking and confinement; two properties observed in Nature. 

Rather than working directly with the quarks and gluons of QCD, another option is to 
construct a model which reproduces the effects of QCD using an effective field theory. This 
is a popular method within the field of nuclear physics, and the route that has be taken for 
this work. More precisely, we utilise a balance between attractive and repulsive meson fields 
to reproduce the binding between fermions that the strong interaction is responsible for. 

1 .2 Neutron Stars 

Although gravity may be the weakest of the four fundamental forces over comparable distance 
scales, it is the most prevalent over (extremely) large distances. It is this force that must be 
overcome for a star to remain stable against collapse. Although a description of this force that 
unifies it with the other three forces has not been (satisfactorily) found. General Relativity 
has proved its worth for making predictions that involve large masses. 

At the time when neutron stars were first proposed by Baade and Zwicky [2], neutrons 
had only been very recently proven to exist by Chadwick [3]. Nonetheless, ever increasingly 
more sophisticated and applicable theories have continually been produced to model the 
interactions that may lead to these incredible structures; likely the most dense configuration 
of particles that can withstand collapse. 

The current lack of experimental data for neutron stars permits a wide variety of mod- 
els [4-9], each of which is able to successfully reproduce the observed properties of neutron 
stars, and most of which are able to reproduce current theoretical and experimental data for 
finite nuclei and heavy-ion collisions [10,11]. The limits placed on models from neutron star 
observations [12-14] do not sufficiently constrain the models, so we have the opportunity to 
enhance the models based on more sophisticated physics, while still retaining the constraints 
above. 

The story of the creation of a neutron star begins with a reasonably massive star, with a 
mass greater than eight solar masses (M > 8 Mq). After millions to billions of years or so 
(depending on the exact properties of the star), this star will have depleted its fuel by fusion 
of hydrogen into ^He, ^He, and larger elements up to iron (the most stable element since has 
the highest binding energy per nucleon). 

At this point, the core of the star will consist of solid iron, as the heaviest elements are 
gravitationally attracted to the core of the star, with successively lighter elements layered 
on top in accordance with the traditional onion analogy. The core is unable to become any 
more stable via fusion reactions and is only held up against gravitational collapse by the 



^At high energies, QCD becomes asymptotically free [,] and can be treated perturbatively. The physics of 
our world however is largely concerned with low energies. 
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1.2. Neutron Stars 



degeneracy pressure of the electrons^. The contents of the upper layers however continue to 
undergo fusion to heavier elements which also sink towards the core, adding to the mass of 
the lower layers and thus increasing the gravitational pressure below. 

This causes the temperature and pressure of the star to increase, which encourages further 
reactions in the upper levels. Iron continues to pile on top of the core until it reaches the 
Chandrasekhar limit of M = 1.4 Mq, at which point the electron degeneracy pressure is 
overcome. The next step is not fully understood"*, but the result is a Type II supernova. 

At the temperatures and pressures involved here, it is energetically favourable for the 
neutrons to undergo /3-decay into protons, electrons (or muons), and antineutrinos according 
to 

n^p+ + e~ + Pe-- (1-1) 

These antineutrinos have a mean-free path of roughly 10 cm [15] at these energies, and are 
therefore trapped inside the star, causing a neutrino pressure bubble with kinetic energy of 
order 10^^ erg = 6.2 x 10^^ MeV [15]. With the core collapsing (and producing even more 
antineutrinos) even the rising pressure of the bubble cannot support the mass of the material 
above and the upper layers begin falling towards the core. 

The sudden collapse causes a shock- wave which is believed to ^ bounce^ at the core and 
expel the outer layers of the star in a mere fraction of a second, resulting in what we know 
as a supernova, and leaving behind the expelled material which, when excited by radiation 
from another star, can be visible from across the galaxy as a supernova remnant (SNR). 

At the very centre of the SNR, the remaining core of the star (naively a sphere of neutrons, 
with some fraction of protons, neutrons and electrons) retains the angular momentum of the 
original star, now with a radius on the order of 10 km rather than 10^ km and thus neutron 
stars are thought to spin very fast, with rotational frequencies of up to 0.716 MHz [IG]. Via a 
mechanism involving the magnetic field of the star, these spinning neutron stars may produce 
a beam of radiation along their magnetic axis, and if that beam happens to point towards 
Earth to the extent that we can detect it, we call the star a pulsar. For the purposes of 
this research, we shall assume the simple case that the objects we are investigating are static 
and non-rotating. Further calculations can be used to extrapolate the results to rotating 
solutions, but we shall not focus on this aspect here. 

A further option exists; if the pressure and temperature (hence energy) of the system 
become great enough, other particles can be formed via weak reactions; for example, hyperons. 
The methods employed in this thesis have the goal of constructing models of matter at super- 
nuclear densities, and from these, models of neutron stars. The outcome of these calculations 
is a set of parameters which describe a neutron star (or an ensemble of them). Of these, 
the mass of a neutron star is an observable quantity. Other parameters, such as radius, 
energy, composition and so forth are unknown, and only detectable via higher-order (or 
proxy) observations. 

The ultimate goal would be finding a physically realistic model based on the interactions 
of particles, such that we are able to deduce the structure and global properties of a neutron 
star based only on an observed mass. This however — as we shall endeavor to show — is easier 
said than done. 



^In accordance with the Pauh Exclusion Principle, no two fermions can share the same quantum state. This 
limits how close two fermions — in this case, electrons — can be squeezed, leading to the degeneracy pressure. 
■^At present, models of supernova production have been unable to completely predict observations. 
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Particle Physics & Quantum Field Theory 



In our considerations of the models that follow we wish to explore ensembles of particles and 
their interactions. In order to describe these particles we rely on Quantum Field Theory 
(QFT), which mathematically describes the 'rules' these particles obey. The particular set of 
rules that are believed to describe particles obeying the strong force at a fundamental level is 
Quantum Chromodynamics (QCD), but as mentioned in the introduction, this construction 
is analytically insolvable, so we rely on a model which simulates the properties that QCD 
predicts. 

In the following sections, we will outline the methods of calculating the properties of 
matter from a field theoretic perspective. 



The first step to calculating any quantity in a Quantum Field Theory is to construct a 
Lagrangian density, which summarizes the dynamics of the system, and from which the 
equations of motion can be calculated. In order to do this, we must define precisely what it 
is that we wish to calculate the properties of. 

The classification schemes of particle physics provide several definitions into which par- 
ticles are identified, however each of these provides an additional piece of information about 
those particles. We wish to describe nucleons (consisting of protons p, and neutrons n) 
which are hadrons^, and are also fermions^. 

We will extend our description to include the hyperons Y (baryons with one or more 
valence strange quarks) consisting of A, S", S*^, S"^, H", and H*^ baryons. The hyperons, 
together with the nucleons, form the octet of baryons (see Fig. B.l). 

We can describe fermions as four-component spinors of plane-wave solutions to the 
Dirac Equation (see later), such that 



where u{p) are four-component Dirac spinors related to plane-waves with wave- vector p that 
carry the spin information for a particle, and which shall be discussed further in Appendix A. 3. 
For convenience, we can group the baryon spinors by isospin group, since this is a degree of 
freedom that will become important. For example, we can collectively describe nucleons as 
a (bi-)spinor containing protons and neutrons, as 



Here we have used the labels for protons and neutrons rather than explicitly using a label for 
isospin. We will further simplify this by dropping the label for spin, and it can be assumed 

^ Bound states of quarks. In particular, bound states of three 'valence' quarks plus any number of quark- 
antiquark pairs (the 'sea' quarks, which are the result of particle anti-particle production via gluons) are called 
baryons. 

^Particles which obey Fermi-Dirac statistics, in which the particle wavefunction is anti-symmetric under 
exchange of particles; the property which leads to the Pauli Exclusion Principle. 



2.1 Lagrangian Density 



ijj = u{p)e 



(2.1) 




(2.2) 
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2. 1. Lagrangian Density 



that this label is implied. We will also require the Dirac Adjoint to describe the antibaryons, 
and this is written as 

ij = V''^7°. (2.3) 

Similarly, we can construct spinors for all the baryons. With these spinors we can con- 
struct a Lagrangian density to describe the dynamics of these particles. Since we are de- 
scribing spin- 1 particles we expect the spinors to be solutions of the Dirac equation which in 
natural units (for which h = c = 1) is written as 

{ij^d^ -M)i; = {i^-M)ip = 0, (2.4) 

and similarly for the antiparticle ^. Feynman slash notation is often used to contract and 
simplify expressions, and is simply defined as^= ^^A^. Here, is the four-derivative, M 
is the mass of the particle, and are the (contravariant) Dirac Matrices, which due to the 
anti-commutation relation of 



2r/"^I, (2.5) 



(where rj = diag(-|-l, —1, —1, —1) is the Minkoswki metric) generate a matrix representation 
of the Clifford Algebra Cl(l, 3). They can be represented in terms of the 2x2 identity matrix 
I, and the Pauli Matrices a, as 

Eq. (2.4) describes free baryons, so we can use this as the starting point for our Lagrangian 
density, and thus if we include each of the isospin groups, we have 

C = Y,4,k{i?>- Mk) A; ke {N, A, S, H}, (2.7) 

k 

where the baryon spinors are separated into isospin groups, as 

V'7v = (^;;), V'A = w, ^P^=(TSj, V'= = (^:';)- (2-8) 

This implies that the mass term is also a diagonal matrix. In many texts this term is simply a 
scalar mass term multiplied by a suitable identity matrix, but that would imply the existence 
of a charge symmetry; that the mass of the proton and of the neutron were degenerate, and 
exchange of charges would have no effect on the Lagrangian density. We shall not make this 
assumption, and will rather work with the physical masses as found in Ref. [17], so Mf. will 
contain distinct values along the diagonal. 

To this point, we have constructed a Lagrangian density for the dynamics of free baryons. 
In order to simulate QCD, we require interactions between baryons and mesons to produce 
the correct phenomenology. Historically, the scalar-isoscalar meson'^ a and vector-isoscalar 
meson iv have been used to this end. Additionally, the vector-isovector p meson has been 
included (for asymmetric matter) to provide a coupling to the isospin channel [18]. 



'Despite it's dubious status as a distinct particle state, rather than a resonance of tttt. 
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In order to describe interactions of the baryons with mesons, we can include terms in the 
Lagrangian density for various classes of mesons by considering the appropriate bilinears that 
each meson couples to. For example, if we wish to include the uj meson, we first observe that 
as a vector meson it will couple to a vector bilinear (to preserve Lorentz invariance) as 



(2.9) 



with coupling strength Qi^, which as we shall see, may be dependent on the baryon that 
the meson is coupled to. The particular coefficients arise from the Feynman rules for meson- 
baryon vertices (refer to Appendix A.l). This particular vertex is written in Feynman diagram 
notation as shown in Fig. 2.1(b). 

This is not the only way we can couple a meson to a baryon. We should also consider the 
Yukawa couplings of mesons to baryons with all possible Lorentz characteristics; for example, 
the u meson can couple to a baryon ?/^, with several different vertices: 



Tp-f^u^ilj, Tpa^i^q^^Tp, and i'q^uj'^ii 



(2.10) 

. The latter two of these 



where represents the baryon four-momentum transfer {qj — qi)^ 
provides a vanishing contribution when considering the mean-field approximation (which shall 
be defined in Section 2.2), since uoo = and q^ = 0, as the system is on average, static. 

If we include the appropriate scalar and vector terms — including an isospin-coupling of 
the p meson — in our basic Lagrangian density we have 



C 



k 



A {it. [id" - Qku^uj" - gpiT(k) -P^)]- Mk + gka<y) i^k ; k e {N, A, S, E}. (2.11) 



The isospin matrices fj^fc) are scaled Pauli matrices of appropriate order for each of the 
isospin groups, the third components of which are given explicitly here as 



r(S)3 







T(A)3 



T-(S)3 















(2.12) 



for which the diagonal elements of T(^k)3 the isospin projections of the corresponding 
baryons within an isospin group defined by Eq. (2.8), i.e. r(p)3 = I^p = 




Fig. 2.1: Interaction vertex for the (a) scalar and (b) vector mesons, where the solid lines 
represent baryons ip, the dashed line represents a scalar meson (e.g. a), and the wavy line 
represents a vector meson (e.g. w^). 
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2.2. Mean-Field Approximation 



At this point it is important that we remind the reader that the conventions in this field 
do not distinguish the use of explicit Einstein summation, and that within a single equation, 
indices may represent summation over several different spaces. To make this clearer, we will 
show an example of a term where all the indices are made explicit; the interaction term for 
the p meson in Eq. (2.11) for which we explicitly state all of the indices 

= ^ 9p 'ipkl'^^k) ■ Pf.ipk 

k 

= EE EEEap (I'U {r^k)Y Pi (^0,' ('-^^^ 

k i,j=l o,/3=l /^=0 a=l 

where here k is summed over isospin groups N, A, S, and H; i and j are summed over flavor 
space (within an isospin group of size fk, e.g. fN = '^,fE = 'i)',Oi and f3 are summed over 
Dirac space; p is summed over Lorentz space; and a is summed over iso-vector space. The 
Pauli matrices, of which (t^)*"' are the elements, are defined in Eq. (2.12). This level of 
disambiguity is overwhelmingly cluttering, so we shall return to the conventions of this field 
and leave the indices as implicit. 

In addition to the interaction terms, we must also include the free terms and field tensors 
for each of the mesons, which are chosen with the intent that applying the Euler-Lagrange 
equations to these terms will produce the correct phenomenology, leading to 

>C = "^ipk [ip [id^ - 9kuj^^ - 9p{T{k) ■ P^)]- Mk + gkacr) ipk 

k 

(2.14) 

where the field tensors for the uj and p mesons are, respectively, 

n^, = d^u:,-d,uj^, R;, = d,,p:-d,pl-gpe-'''plp:. (2.15) 

This is the Lagrangian density that we will begin with for the models we shall explore herein. 
Many texts (for example, Refs. [19-22]) include higher-order terms {0{a^), O(cr^), . . .) and 
have shown that these do indeed have an effect on the state variables, but in the context of 
this work, we shall continue to work at this order for simplicity. It should be noted that the 
higher order terms for the scalar meson can be included in such a way as to trivially reproduce 
a framework consistent with the Quark-Meson Coupling model that shall be described later, 
and thus we are not entirely excluding this contribution. 

2.2 Mean- Field Approximation 

To calculate properties of matter, we will use an approximation to simplify the quantities 
we need to evaluate. This approximation, known as a Mean-Field Approximation (MFA) is 
made on the basis that we can separate the expression for a meson field a into two parts: a 
constant classical component, and a component due to quantum fluctuations; 



— Q^classical ~l~ Ctquantum- 



(2.16) 



2.3. Symmetries 
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If we then take the vacuum expectation value (the average value in the vacuum) of these 
components, the quantum fluctuation term vanishes, and we are left with the classical com- 
ponent 

(a) = (Oclassical)- (2-17) 

This component is what we shall use as the meson contribution, and we will assume that 
this contribution (at any given density) is constant. This can be thought of as a background 
'field' on top of which we place the baryon components. For this reason, we consider the case 
of infinite matter, in which there are no boundaries to the system. The core of lead nuclei 
(composed of over 200 nucleons) can be thought of in this fashion, since the effects of the 
outermost nucleons are minimal compared to the short-range strong nuclear force. 

Furthermore, given that the ground-state of matter will contain some proportion of proton 
and neutron densities, any flavor-changing meson interactions will provide no contribution in 
the MFA, since the overlap operator between the ground-state \^) and any other state |^) is 
orthogonal, and thus 

( ^ I O = ^n- (2-18) 

For this reason, any meson interactions which, say, interact with a proton to form a neutron 
will produce a state which is not the ground state, and thus provides no contribution to 
the MFA. We will show in the next section that this is consistent with maintaining isospin 
symmetry. 

2.3 Symmetries 

In the calculations that will follow, there are several terms that we will exclude from our 
considerations ab initio (including for example, some that appear in Eq. (2.14)) because they 
merely provide a vanishing contribution, such as the quantum fluctuations mentioned above. 
These quantities shall be noted here, along with a brief argument supporting their absence 
in further calculations. 

2.3.1 Rotational Symmetry and Isospin 

The first example is simple enough; we assume rotational invariance of the fields to conserve 
Lorentz invariance. In order to maintain rotational invariance in all frames, we require that 
the spatial components of vector quantities vanish, leaving only temporal components. For 
example, in the MFA the vector-isoscalar meson four- vector uj^ can be reduced to the temporal 
component ojq, and for notational simplicity, we will often drop the subscript and use (a) for 
the a meson mean-field contribution. 

A corollary of the MFA is that the field tensor for the rho meson vanishes; 

R^u = d,pl - d^pl - gpe'^'^p'lpl —> i?So = 0, (2.19) 

Mr A 

since the derivatives of the constant terms vanish and (pb x Po) = 0. The same occurs for the 
omega meson field tensor 



^flU — dnUJy — dpUln > — 0. 

^ MFA 



(2.20) 
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2.3. Symmetries 



We also require rotational invariance in isospin space along a quantization direction of z = 3 
(isospin invariance) as this is a symmetry of the strong interaction, thus only the neutral 
components of an isovector have a non-zero contribution. This can be seen if we examine the 
general 2x2 unitary isospin transformation, and the Taylor expansion of this term 

V'(x) ^ i:\x) = e'^-^/'^i^{x) (l + if- 6/2] ^^{x), (2.21) 

|e|<i V / 

where ^ = (^i, ^3) is a triplet of real constants representing the (small) angles to be rotated 
through, and r are the usual Pauli matrices as defined in Eq. (2.12). As for the p mesons, 
we can express the triplet as linear combinations of the charged states, as 

P = {pi,P2,p3) = (^-^{p+ + P-), -^(p_ -/)+), po^ . (2.22) 
The transformation of this triplet is then 

p{x) ^ p\x) = e'^-^p{x). (2.23) 

where {T^)jk = —i^ijk is the adjoint representation of the SU(2) generators; the spin-1 Pauli 
matrices in the isospin basis, a.k.a. S0(3). We can perform a Taylor expansion about 6 = 0, 
and we obtain 

Pj{x) [djk + i{T')jkOi\pk{x). (2.24) 

We can therefore write the transformation as 

p(x) — > p(x)-6xp(x). (2.25) 

Writing this out explicitly for the three isospin states, we obtain the individual transformation 
relations 

p (x) ^ p (x) - 6 X p (x) = {pi - 62P3 + 63P2, P2 - 6iP3 + 63P1, p3 - 61P2 + 62P1). 

(2.26) 

If we now consider the rotation in only the z = 3 direction, we see that the only invariant 
component is p^ 

p (x) > (pi + 63P2, P2 + 93Pi, Ps) ■ (2.27) 

91=0 
92=0 

If we performed this rotation along another direction — i.e. 1, 2, or a linear combination 
of directions — we would find that the invariant component is still a linear combination of 
charged states. By enforcing isospin invariance, we can see that the only surviving p meson 
state will be the charge-neutral state pa = po- 



2.3.2 Parity Symmetry 

We can further exclude entire isospin classes of mesons from contributing since the ground- 
state of nuclear matter (containing equal numbers of up and down spins) is a parity eigenstate, 
and thus the parity operator V acting on the ground-state produces 



V\0) = ±\0). 



(2.28) 



2.3. Symmetries 
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Noting that the parity operator is idempotent (7^^ = I), inserting the unity operator into the 
ground-state overlap should produce no effect; 

{o\o) = {o\i\o) = {o\rr\o) = {o\{±f\o) = {o\o). (2.29) 

We now turn our attention to the parity transformations for various bilinear combinations 
that will accompany meson interactions. For Dirac spinors ip{x) and 4^{x) the parity trans- 
formation produces 

V^p{t,x)V = -f'^iP{t,-x), ^^^^^ 
'Ptp{t,x)V = 4){t,-x)-f^, 

where we have removed the overall phase factor exp since this is unobservable and can 
be set to unity without loss of generalisation. We can also observe the effect of the parity 
transformation on the various Dirac field bilinears that may appear in the Lagrangian density. 
The five possible Dirac bilinears are: 

^/;^/;, V'T'^V', ^V'[7^,7l'0, V't'^T^V', i^J^ip, (2.31) 

for scalar, vector, tensor, pseudo-vector and pseudo-scalar meson interactions respectively, 
where 75 is defined as 

75 = «7o7i7273 = ( n ) ' (^•"^^) 

in the commonly used Dirac basis. By acting the above transformation on these bilinears we 
obtain a result proportional to the spatially reversed wavefunction ?/^(t, —x), 

Vi^^l^V = +i^i){t,-x), (2.33) 

' ^' \ -i^-l^^ijit^-x) for /i = 1,2, 3, ^ > 

(^^,5,1,-0 - / -V'7VV'(i, ) for/i = 0. 



' ' [ +'>p-f>^-f^'>p{t,-x) for /i = 1,2, 3, ^ ' 

Vii/j-f^tpV = -itp-f^i;{t,-x). (2.36) 

By inserting the above pseudo-scalar and pseudo- vector bilinears into the ground-state overlap 
as above, and performing the parity operation, we obtain a result equal to its negative, and 
so the overall expression must vanish. For example 

{0\iij-f^ip\0) = {0\Vii)-i^il)V\0) = {0\ - ii)-i^i^\0) = 0. (2.37) 

Thus all pseudo-scalar and pseudo- vector meson contributions — such as those corresponding 
to TT and K — provide no contribution to the ground-state in the lowest order. We will show 
later in Chapter 3.5 that mesons can provide higher order contributions, and the pseudo- 
scalar TT mesons are able to provide a non-zero contribution via Fock terms, though we will 
not calculate these contributions here. 
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2.4. Fermi Momentum 



2.4 Fermi Momentum 

Since we are dealing with fermions that obey the Pauh Exclusion Principle'*, and thus Fermi- 
Dirac statistics^, there will be restrictions on the quantum numbers that these fermions may 
possess. When considering large numbers of a single type of fermion, they will each require a 
unique three-dimensional momentum k since no two fermions may share the same quantum 
numbers. 

For an ensemble of fermions we produce a 'Fermi sea' of particles; a tower of momentum 
states from zero up to some value 'at the top of the Fermi sea'. This value — the Fermi 
momentum — will be of considerable use to us, thus it is denoted kp- 

Although the total baryon density is a useful control parameter, many of the parameters of 
the models we wish to calculate are dependent on the density via kp- The relation between the 
Fermi momentum and the total density is found by counting the number of momentum states 
in a spherical volume up to momentum kp (here, this counting is performed in momentum 
space). The total baryon density — a number density in units of baryons/fm'^ , usually denoted 
as just fm~'^ — is simply the sum of contributions from individual baryons, as 

Ptotai = Ep^ = E ^^70^ I ^(^^» -\k\)d'k = Y,^, (2.38) 

i i i 

where here, i is the set of baryons in the model, Jj is the spin of baryon i (where for the 
leptons and the octet of baryons, Jj = ^), and 9 is the Heaviside step function defined as 

which restricts the counting of momentum states to those between and kp. 

We define the species fraction for a baryon lepton £, or quark q as the density fraction 
of that particle, denoted by 5^, such that 

Y, = -^- ie{B,^,q}. (2.40) 

Ptotal 

Using this quantity we can investigate the relative proportions of particles at a given total 
density. 

2.5 Chemical Potential 

In order to make use of statistical mechanics we must define the some important quantities. 
One of these will be the chemical potential ;U, also known as the Fermi energy ep] the energy 
of a particle at the top of the Fermi sea, as described in Appendix A. 5. 5. This energy is the 
relativistic energy of such a particle, and is the energy associated with a Dirac equation for 
that particle. For the simple case of a non-interacting particle, this is 

I^B = ep^ = ^kl^+Ml. (2.41) 
^That no two fermions can share a single quantum state. 

®The statistics of indistinguishable particles with half-integer spin. Refer to Appendix A. 5. 5. 



2.5. Chemical Potential 
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In the case that the baryons are involved in interactions with mesons, we need to introduce 
scalar and temporal self-energy terms, which (for example) for Hartree- level QHD using a 
mean-field approximation are given by 

^% = -9Ba{a), ^% = 9bM + 9phB{p), (2.42) 

where I^b is the isospin projection of baryon B, defined by the diagonal elements of Eq. (2.12), 
and where the scalar self-energy is used to define the baryon effective mass as 

M*B = Mb + S|j = Mb - qb. {ct) , (2.43) 

These self-energy terms affect the energy of a Dirac equation, and thus alter the chemical 
potential, according to 

^iB = ^kl^ + {MB + T.%Y + . (2.44) 

Eq (2.43) and Eq. (2.44) define the important in-medium quantities, and the definition of 
each will become dependent on which model we are using. 

For a relativistic system such as that which will consider here, each conserved quantity 
is associated with a chemical potential, and we can use the combination of these associated 
chemical potentials to obtain relations between chemical potentials for individual species. In 
our case, we will consider two conserved quantities: total baryon number and total charge, 
and so we have a chemical potential related to each of these. We can construct the chemical 
potential for each particle species by multiplying each conserved charge by its associated 
chemical potential to obtain a general relation. Thus 

Hi = Bifin - Qi^i-e, (2.45) 

where; i is the particle species (which can be any of the baryons) for which we are constructing 
the chemical potential; Bi and Qi are the baryon number ('baryon charge', which is unitless) 
and electric charge (normalized to the proton charge) respectively; and /i„ and /ie are the 
chemical potentials of neutrons and electrons, respectively. Leptons have i?£ = 0, and all 
baryons have Bb = +1. The relations between the chemical potentials for the octet of 
baryons are therefore derived to be 

/"A = /^SO = /iHO = /in, 

/^E- = A'H- = /in + /^e, ^2 46) 

Mp = A'E+ = — ^J'e, 

IJ^tJi = fJ-e- 

A simple example of this is to construct the chemical potential for the proton (for which the 
associated charges are Bp = +1 and Qp = +1); 

/Up = /^n - /^e- (2.47) 
This can be rearranged to a form that resembles neutron /?-decay 

fin = fJ'p + fJ'e- (2.48) 



If we were to consider further conserved charges, such as lepton number for example, 
we would require a further associated chemical potential. In that example, the additional 
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2.6. Explicit Chiral Symmetry (Breaking) 



chemical potential would be for (anti) neutrinos /ip. The antineutrino would be required to 
preserve the lepton number on both sides of the equation; the goal of such an addition. Since 
we shall consider that neutrinos are able to leave the system considered, we can ignore this 
contribution ab inito. The removal of this assumption would alter Eq. (2.48) to include the 
antineutrino, as would normally be expected in /3-decay equations 

IJ-n = IJ'p + IJ-e + fJ'U- (2.49) 



2.6 Explicit Chiral Symmetry (Breaking) 

One of the most interesting symmetries of QCD is chiral symmetry. If we consider the QCD 
Lagrangian density to be the sum of quark and gluon contributions, then in the massless 
quark limit {mg = 0); 

■^QCD = Cg + Cq 

= -^G;i,Gr + V5.i7^a^V'.-5^>O^T;^V', , (2.50) 

where here, ilJi{x) is a quark field of color i G {r, 5, 6}, A'^{x) is a gluon field with color 
index a S {1, . . . , 8}, Tj" is a generator*^ for SU(3), g is the QCD coupling constant, and GJJ^ 
represents the gauge-invariant gluonic field strength tensor, given by 

Gl, = [d,,Al]-gr^-AlAl, (2.51) 

written with the structure constants Left- and right-handed components of Dirac fields 
can be separated using the projection operators 

= i^V', (2.52) 

using the definition of 75 of Eq. (2.32), and so the quark terms in the QCD Lagrangian 
density (the gluon terms are not projected) can be written in terms of these components as 

CSp = i^^D^-i^^f + i^^^D^^^4i^ . (2.53) 

This Lagrangian density is invariant under rotations in U(l) of the left- and right-handed 
fields 

U(1)l: V'L^e^^^V'L, ^R^^R, (2.54) 
U(1)r: ^R^e^"«VR, V'L^^L, (2.55) 

where and or are arbitrary phases. This invariance is the chiral \J (1) l <S>V (1) r symmetry. 
The Noether currents associated with this invariance are then 

= i^Ll'^i^L, 4 = ^Fr^^I^r, (2.56) 



^For example, T" = X"^ /2 using the Hermitian Gell-Mann matrices \a 
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and as expected, these currents are conserved, such that 5^ Jj^ = dfj,J^ = according to the 
Dirac Equation. These conserved currents can be alternatively written in terms of conserved 
vector and axial- vector currents, as 

J'l = J'r = (2.5T) 

where here, V^^ and denote the vector and axial-vector currents respectively — the dis- 
tinction of A^^ here from the gluon fields in Eq. (2.50) is neccesary — and these are defined 
by 

V^' = iPjf'i;, Af" = V't'^TsV', (2-58) 

and which are also conserved, thus dfj,V^ = d^A^ = 0. The chiral symmetry of 11(1)1, (8) U(l)ij 
is therefore equivalent to invariance under transformations under U(l)v/ ® U(l)yi, where we 
use the transformations 

U(l)y : V ^ e^^^V', V' ^ V''^e-^°^7o, (2.59) 
U(1)a : ^ ^ e*"-^^5^, ^ V'^e-^"-^T5 7o. (2.60) 

Using the anticommutation relation 

{75, 7m} = 757m + 7m75 = (2.61) 

we can evaluate the effect that the vector and axial- vector transformations have on the QCD 
Lagrangian density, and we find that both transformations are conserved. If we now consider 
a quark mass term Cm in the QCD Lagrangian density, the fermionic part becomes 

^QCD = ^q + ^m = i^i {i-i^{D^)ij - mSij) (2.62) 

For the purposes of these discussions, we can set the masses of the quarks to be equal without 
loss of generality. Although the massless Lagrangian density possesses both of the above 
symmetries, the axial vector symmetry — and hence chiral symmetry — is explicitly broken 
by this quark mass term; 

jCm = —'iprmp ^-^^^ —^TTji^^'^o'A^ —tpmip. (2.63) 
The vector symmetry is nonetheless preserved when including this term. 



2.7 Dynamical Chiral Symmetry (Breaking) 

Even with a massless Lagrangian density, it is possible that chiral symmetry becomes dy- 
namically broken, and we refer to this as Dynamically Broken Chiral Symmetry, or DCSB. 

Following the description of Ref. [23], if we consider the basic Lagrangian density of QCD 
to be 

Cqcb = (iriD^y^j - m5,,) V, - ^G^.G^, (2.64) 
with the definitions as in the previous section, of 

G^, = 9^^^-9,^^-<7r^^4^^, D, = d, + igA'^T-, (2.65) 
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with standard definitions of other terms, then we can write the sum of all QCD One-Particle 
Irreducible (1-PI) diagrams^ with two external legs as shown in Fig. 2.2; illustrating the 
quark self-energy. The expression for the renormalized quark self-energy in d dimensions is 

-iS(p) = prg'l -0^,{il,){iS{qmD^''{v-q)WM,P)), (2.66) 

where Z^. is a renormalization constant, g is the QCD coupling, and q is the loop momentum. 

In the absence of matter fields or background fields (the Lorentz-covariant case), we can 
write this self-energy as a sum of Dirac- vector and Dirac-scalar components, as 

S(p) = ^Edv(/) + Sds(p')- (2.67) 

where Edv(j'^) is the Dirac-vector component, and Sds(?'^) is the Dirac-scalar component. 
These must both be functions of p^, since there are no other Dirac-fields to contract with, 
and Ti{p) is a Lorentz invariant quantity in this case. 

For the purposes of our discussion in this section, we will approximate the Dirac-vector 
component of the self-energy to be Sdv ~ 1, in which case the self-energy is dependent only 
on the Dirac-scalar component. 

Even with a massless theory (m = 0) it is possible that the renormalized self-energy 
develops a non-zero Dirac-scalar component, thus Sds(p^) 7^ 0. This leads to a non-zero 
value for the quark condensate {ipqipq), and in the limit of exact chiral symmetry, leads to 
the pion becoming a massless Goldstone boson. Thus chiral symmetry can be dynamically 
broken. With the addition of a Dirac-scalar component of the self-energy, the Lagrangian 
density becomes 

JI^QCD = {ir{Df,)ij - (m + Eds)<5.j) - ^G^^Cr, (2.68) 



iD 




Fig. 2.2: Feynman diagram for the QCD self-energy for a quark, as given by the Dyson- 
Schwinger Equation (DSE). The full expression for this is given in Eq. (2.66). 



Diagrams that cannot be made into two separate disconnected diagrams by cutting an internal line are 
called One-Particle Irreducible, or 1-PI. 



2.8. Equation of State 
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and we can define a dynamic quark mass via the gap equation; 

m* = m + Eds- (2.69) 

We will continue this discussion in Section 3.4, in which we will describe a particular model 
for Eds in order to describe DCSB. 



2.8 Equation of State 



In order to investigate models of dense matter, we need to construct an Equation of State 
(EOS), which is simply a relation between two or more state variables — those which thermo- 
dynamically describe the current state of the system, such as temperature, pressure, volume, 
or internal energy — under a given set of physical conditions. With this, we will be able to 
investigate various aspects of a model and compare differences between models in a consistent 
fashion. 

For our purposes, we use the total baryon density /Ototai as the control parameter of this 
system, and so we need to obtain the connection between, say, the energy density £, the 
pressure P, and this total baryon density, i.e. 

= ^^(Ptotal), P = P(ptotal). (2.70) 

State variables are important quantities to consider. Within any transition between states 
the total change in any state variable will remain constant regardless of the path taken, since 
the change is an exact differential, by definition. For the hadronic models described herein, 
the EOS are exact, in that they have an analytic form; 

2 d f <f (ptotal) 



^'(Ptotai) = Ptot.ij^ • (2.71) 

Cptotal V Ptotal / 

As simple as this exact form may seem, the derivative complicates things, and we will find 
it easier to calculate the pressure independently. Nonetheless, this expression will hold true. 
More interestingly, this expression is equivalent to the first law of thermodynamics in the 
absence of heat transfer; i.e. 

PdV = -dE, (2.72) 

(the proof of which can be found in Appendix A. 5. 6) which assures us that the theory is 
thermodynamically consistent. 

A notable feature of each symmetric matter EOS we calculate is the effect of saturation; 
whereby the energy per baryon for the system possesses a global minimum at a particular 
value of the Fermi momentum. This can be considered as a binding energy of the system. 
In symmetric matter (in which the densities of protons and neutrons are equal), the nucleon 
Fermi momenta are related via kp = kp^ = kpp^ and the energy per baryon (binding energy) 
E is determined via 

E= -^(e-y^pBMB] . (2.73) 

Ptotal \ B J 

In order to reproduce (a chosen set of) experimental results, this value should be an ex- 
tremum of the curve with a value of Eq = —15.86 MeV at a density of po = 0.16 fm~^ (or 
the corresponding Fermi momentum kp^^). 
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2.9. Phase Transitions 



The nucleon symmetry energy agym is approximately a measure of the energy difference 
between the energy per baryon (binding energy) of a neutron-only model and a symmet- 
ric nuclear model (essentially a measure of the breaking of isospin symmetry). A more 
formal expression (without assuming degeneracy between nucleon masses, as derived in Ap- 
pendix A. 5. 4) is 



^sym 



37r2 



r4 + 
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k'p 
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kl + (M*)2 



(2.74) 



At saturation, this should take the value of (asym)o = 32.5 MeV (for an analysis of values, 
see Ref. [24]). 

Another important aspect of an EOS is the compression modulus K which represents the 
stiffness of the EOS; the ability to withstand compression. This ability is intimately linked 
to the Pauli Exclusion Principle in that all other things being equal, a system with more 
available states (say, distinguishable momentum states) will have a softer EOS, and thus a 
smaller compression modulus. The compression modulus itself is defined as the curvature of 
the binding energy at saturation, the expression for which is 
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Ptotal 



9 



2 

Ptotal 



d2 



dplotaX VPtotal, 



(2.75) 



P=PO 



The motivation for this is that by compressing the system, the energy per baryon will rise. 
The curvature at saturation determines how fast that rise will occur, and thus how resistant 
to compression the system is. Experimentally, this is linked to the properties of finite nuclei, 
particularly those with a large number of nucleons, and the binding of these within a nucleus. 

According to Ref. [18] this should have a value in the range 200-300 MeV, and we will 
calculate the value of K for each of the models to follow for comparison. 



2.9 Phase Transitions 

In order to consider transitions between different phases of matter we must use statistical 
mechanics. The simplest method of constructing a phase transition — known as a 'Maxwell 
transition' — is an isobaric (constant pressure) transition constructed over a finite density 
range. A transition of this form remains useful in understanding the liquid-gas style phase 
transition that occurs within QHD, which is a first-order transition (similar to that of ice 
melting in a fiuid) with the phases being separated by a non-physical negative-pressure region. 
The inclusion of a Maxwell transition to this simple model for QHD removes this unphysical 
region and replaces it with a constant pressure phase. 

The method for constructing a Maxwell transition will not be covered here, though in- 
depth details can be found in Ref. [25]. We can however extract the transition densities from 
Ref. [18] to reproduce the results, which are shown later in Fig. 5.5 for the various varieties 
of QHD. The more sophisticated method of constructing a phase transition — the 'Gibbs 
transition' [26] that we have used for the results produced herein — relies on a little more 
statistical mechanics. A comparison between the Maxwell and Gibbs methods for models 
similar to those used in this work can be found in Ref. [27]. For a full in-depth discussion of 
this topic, see Ref. [28]. 
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If we consider a homogeneous (suitable for these mean-field calculations) system with en- 
ergy E, volume V, Nm particles of type m, and entropy S which depends on these parameters 
such that 

S = S{E,V,Ni,...,Nm), (2.76) 

then we can consider the variation of the entropy in the system as a function of these param- 
eters, resulting in 

J V,Ni,...,N^ / E,Ni,...,Nm i=i \^^^^/V,Nj^i 

using standard statistical mechanics notation whereby a subscript X on a partial derivative 
(dA/dB)-^ denotes that X is explicitly held constant. Eq. (2.77) should be equal to the 
fundamental thermodynamic relation when the number of particles is fixed, namely 

,S = ^ = f?™:. (2.T8) 

Here, the symbol d denotes the inexact differential, since the heat Q is not a state function 
— does not have initial and final values — and thus the integral of this expression is only true 
for infinitesimal values, and not for finite values. Continuing to keep the number of each type 
of particle Ni constant, a comparison of coefficients between Eqs. (2.77) and (2.78) results in 
the following relations: 

dS\ 1 /dS\ P , ^ , 

(2.79) 



J V,Ni,-,Nm ^ \(^^ J E,N„...,N^ ^ 

To provide a relation similar to Eq. (2.79) for the case where dNi ^ 0, one defines fij — the 
chemical potential per molecule — as 

We can now re-write Eq. (2.77) with the definitions in Eq. (2.79) for the case where the 
particle number can change, as 

P ™ 

dS = -dE + -dV-Y, Y'^N,, (2.81) 
1=1 

which can be equivalently written in the form of the fundamental thermodynamic relation 
for non-constant particle number, 

m 

dE = TdS - PdV + ^ ^l^dNi. (2.82) 

i=l 

If we now consider a system X of two phases A and then we can construct relations 
between their parameters by considering the following relations: 

Ex = Ea + Eb, 

Vx = Va + Vb, (2.83) 
Nx = Na + Nb- 
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If we consider that these quantities are conserved between phases, we find the following 
conservation conditions 

dEx = ^ dEA + dEs = ^ dEA = -dEs, 

dVx = dVA + dVe = ^ dVA = -dVu, (2.84) 

dNx = ^ dNA + dNB = ^ dNA = -dNs- 

The condition for phase equilibrium for the most probable situation is that the entropy 
must be a maximum for 5 = SxiEx, Vx, A^x) = S{Ea, Va, Na; Eb, Vb, Nb), which leads to 

dSx = dSA + dSB = 0. (2.85) 

Thus, inserting Eq. (2.81) we find 

dS = l^-^dEA + ^dVA - ^dNA^ + {-f^dEB + ^dVB - diVs) • (2.86) 

If we now apply the result of Eq. (2.84) we can simplify this relation to 

dS = 0=(^-^)dEA+(^-^)dVA-(!;i^-!^)dNA (2.87) 
\Ta TbJ \Ta TbJ \Ta TbJ ^ ^ 

and thus for arbitrary variations of Ea, Va and Na, each bracketed term must vanish sepa- 
rately, so that 

J_ ^ J_ Pa ^ Pb_ ^lA ^ P;B_ 

-LA -LB -LA -LB -LA -L B 

Eq. (2.88) implies that at the phase transition the system will be isentropic {dS = 0), isother- 
mal {dT = 0), isobaric {dP = 0), and isochemical {dfi = 0), where the terms S, T, V, and fi 
now refer to the mean values rather than for individual particles. 

We only require two systems at any one time when considering a mixture of phases, for 
example, a neutron ('neutron phase') can transition to a proton and an electron ('proton and 
electron phase') provided that the condition //n = /^p + fJ-e is met. 

For a phase transition between hadronic- and quark-matter phases then, the conditions 
for stability are therefore that chemical, thermal, and mechanical equilibrium between the 
hadronic H, and quark Q phases is achieved, and thus that the independent quantities in 
each phase are separately equal. Thus the two independent chemical potentials (as described 
in Sec. 2.5) /x„ and /ig are each separately equal to their counterparts in the other phase, 
i.e. [{Hn)H = {fJ'n)Q], and [{f-ie)H = {fJ'e)Q] for chemical equilibrium; [Th = Tq] for thermal 
equilibrium; and [Ph = Pq] for mechanical equilibrium. 

An illustrative example of these relations is shown in Fig. 2.3 in which the values of the 
independent chemical potentials and fie, as well as the pressure P for a hadronic phase and 
a quark phase are plotted for increasing values of total density ptotai- In this case, the quark 
matter data is calculated based on the hadronic matter data, using the chemical potentials 
in the hadronic phase as inputs for the quark phase calculations, and as such the chemical 
potentials are — by construction — equal between the phases. As this is an illustrative example 
of the relations between the phases, no constraints have been imposed to reproduce a phase 
transition yet. 

In this figure, the low-density points correspond to small values of fin, and we see that 
for densities lower than some phase transition density ptotai < Ppt the hadronic pressure is 
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Fig. 2.3: (Color Online) Illustrative locus of values for the independent chemical potentials 
/ig and fin, as well as the pressure P for phases of hadronic matter and deconfined quark 
matter. Note that pressure in each phase increases with density, and that a projection onto 
the /in/^e plane is a single line, as ensured by the chemical equilibrium condition. 



greater than the quark pressure and thus the hadronic phase is dominant. At the transition 
the pressures are equal, and thus both phases can be present in a mixed phase, and beyond 
the transition the quark pressure is greater than the hadronic pressure indicating that the 
quark phase becomes dominant. 

Note that for all values of the total density, the chemical potentials in each phase are 
equal, as shown by the projection onto the /i^/^e plane. 

In our calculations, we will only investigate these two phases independently up to the 
phase transition, at which point we will consider a mixed phase, as shall be described in the 
next section. 

We consider both phases to be cold on the nuclear scale, and assume T = so the 
temperatures are also equal, again by construction. We must therefore find the point — if it 
exists — at which, for a given pair of independent chemical potentials, the pressures in both 
the hadronic phase and the quark phase are equal. 

To find the partial pressure of any baryon, quark, or lepton species i we use 

P^=(l^±4^f^^lM^d% (2.89) 



3(27r)3 



+ [M* 



where the number of colors is Nc = 3 for quarks, Nc = 1 for baryons and leptons, and where 9 
is the Heaviside step function defined in Eq. (2.39). To find the total pressure in each phase. 
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we sum the pressures contributions in that phase. The pressure in the hadronic phase H is 
given by 

PH = Y.^j + 12^^ + Yl (2.90) 

j £ m 

in which j represents the baryons, i represents the leptons, m represents the mesons appearing 
in the particular model being considered, if they appear, and the pressure in the quark phase 
Q is given by 

^Q = E^'' + E^^-^' (2.91) 

where q represents the quarks, and B denotes the bag energy density, which we shall discuss 
further in Section 3.3. 

In order to determine the EOS beyond the point at which the pressures are equal, we 
need to consider the properties of a mixed phase. 

2.10 Mixed Phase 

Once we have defined the requirements for a phase transition between two phases, we must 
consider the possibility of a mixed phase (MP) containing proportions of the two phases. 
This adds a further degree of sophistication to a model; we can not only find equations of 
state for hadronic matter and quark matter and simply stitch them together, but we can also 
allow the transition between these to occur gradually. 

To calculate the mixed phase EOS, we calculate the hadronic EOS with control parameter 
Ptotai ) and use the independent chemical potentials /i„ and He as inputs to determine the quark 
matter EOS, since we can determine all other Fermi momenta given these two quantities. We 
increase /Ototai until we find a density — if it exists — at which the pressure in the quark phase 
is equal the pressure in the hadronic phase (if such a density cannot be found, then the 
transition is not possible for the given models). 

Assuming that such a transition is possible, once we have the density and pressure at 
which the phase transition occurs, we change the control parameter to the quark fraction x 
(which is an order parameter parameterizing the transition to the quark matter phase) which 
determines the proportions of hadronic matter and quark matter. If we consider the mixed 
phase to be composed of some fraction of hadronic matter and some fraction of quark matter, 
then the mixed phase of matter will have the following properties: the total density will be 

PMP = (1-x) PHP + XPQP, (2.92) 

where php and /9qp are the densities in the hadronic and quark phases, respectively. A factor 
of three in the equivalent baryon density in the quark phase, 

PCIP = \'^Pq = {pu + Pd + Ps)/'^, (2.93) 

arises because of the restriction that a baryon contains three quarks. 

According to the condition of mechanical equilibrium detailed earlier, the pressure in the 
mixed phase will be 

Pmp = Pm> = Pqp. (2.94) 



2. 1 1. Stellar Matter 
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We can step through values < x < 1 and determine the properties of the mixed phase, 
keeping the mechanical stability conditions as they were above. In the mixed phase we need 
to alter our definition of charge neutrality; while previously we have used the condition that 
two phases were independently charge-neutral, such as n — )• p+ + e~ , it now becomes possible 
that one phase is (locally) charged, while the other phase carries the opposite charge, making 
the system globally charge-neutral. This is achieved by enforcing 

= a-x) PRP + XPQP + Pi, (2.95) 

where this time we are considering charge-densities, which are simply the sum of densities 
multiplying their respective charges 

Pt = Y.Q^Pr^ zG{HP,QP,^}, (2.96) 

j 

where j are the all individual particles modelled within the grouping i. For example, the 

quark charge-density in a non-interacting quark phase is given by 

^ 2 11 

^QP = Yl ^'iP'i = 3^^" ~ '^P'^ ~ 3^*- 
<? 

We continue to calculate the properties of the mixed phase for increasing values of x 
until we reach x = 1) at which point the mixed phase is now entirely charge-neutral quark 
matter. This corresponds to the density at which the mixed phase ends, and a pure quark 
phase begins. We can therefore continue to calculate the EOS for pure charge-neutral quark 
matter, once again using /Ototai as the control parameter, but now where the total density is 
the equivalent density as defined in Eq. (2.93). 



2.11 Stellar Matter 

The equations of state described above are derived for homogeneous infinite matter. If we 
wish to apply this to a finite system we must investigate the manner in which large ensembles 
of particles are held together. The focus of this work is 'neutron stars', and we must find a 
way to utilise our knowledge of infinite matter to provide insight to macroscopic objects. For 
this reason, we turn to the theory of large masses; General Relativity. 

The Tolman-Oppenheimer-Volkoff (TOV) equation [29] describes the conditions of sta- 
bility against gravitational collapse for an EOS, i.e. in which the pressure gradient is sufficient 
to prevent gravitational collapse of the matter. The equations therefore relate the change 
in pressure with radius to various state variables from the EOS. To preserve continuity, the 
equations are solved under the condition that the pressure at the surface of the star must be 
zero. 

The TOV equation is given by 

dP _ G{P/^ + £){M[r) + Ar^^P/f)_ ^^^^^ 



or, in Planck units^ 



dr r{r -2GM{r)/c'^) 

dP _ {P + £) (M(r) + Wp) 
'dr ~ r(r - 2M(r)) 



(2.99) 



'in which certain fundamental physical constants are normalized to unity, viz h = c = G = 1. 
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where the mass within a radius R is given by integrating the energy density, as 



M{R) = [ 



R 



47rr^i5(r) dr. 



(2.100) 



Jo 



For a full derivation of these equations, refer to Appendix A. 4. 

Supplied with these equations and a derived EOS, we can calculate values for the total 
mass and total radius of a star^ for a given central density. We refer to these values as the 
'stellar solutions.' 

This is particularly interesting, since the mass of a neutron star is observable (either 
via observing a pair of objects rotating about a barycenter^*^, or some other indirect/proxy 
measurement), yet the radius is not directly observable, as stars are sufficiently distant that 
they all appear as 'point-sources'. With these calculations, we produce a relationship between 
two quantities: the stellar mass and the stellar radius, of which only the mass is currently 
observable, and even this is not always so. This provides useful data for further theoretical 
work requiring both quantities, as well as an opportunity to place theoretical bounds on 
future experimental observations. 

In addition to this data, since we are able to solve our equations for the radial distance 
from the centre of the star, we can provide data that current experiments can not; we can 
investigate the interior of a neutron star, by calculating the proportions of various particles 
at successive values of internal radius and/or density. This allows us to construct a cross- 
section of a neutron star, investigate the possible contents, and examine the effects that 
various changes to the models have on both the internal and external properties. 



We have noted earlier that the coupling of baryons to mesons is dependent on isospin group. 
The physics leading to this result is highly non-trivial, but is often neglected in the literature. 
We will therefore outline the process involved in determining the relations between baryon- 
meson couplings. 

In order to determine the normalized relations between the point vertex couplings of 
various mesons to the full baryon octet gBm it is common to express the octet as a 3 x 3 
matrix in flavor space as 



This has been constructed as an array where rows and columns are distinguished by rotations 



^Stellar objects in these calculations are assumed to be static, spherically symmetric, and non-rotating, 
as per the derivation of this equation. For studies of the effect of rapid rotation in General Relativity see 



A common centre of mass for the system, the point about which both objects will orbit, which is the 
balance point of the gravitational force. In this case, the mass measurements are simplified. 



2.12 SU(6) Spin-Flavor Baryon-Meson Couplings 




(2.101) 



Refs. [30,31]. 
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in flavor space, which can be seen if we observe the quark content of these baryons; 

u uus s — )• d uud \ 



B 



I uds d 
u^d 

dds 
d — s 

\ ssd 



dus 



uss 



ddu 
sud ) 



(2.102) 



The vector meson octet (J^ = 1~) can be written in a similar fashion as 



pvec 
oct 



K*- 



p+ K*+ 



\ 



(2.103) 



^6 j 



which, along with the singlet state, P^^^^ = -^diag(a;o, wq, wq) defines the vector meson nonet 

pvec ^ pvec ^pvcc^ (2.104) 

Furthermore, the scalar meson octet (J^ = 0+) can be written as 



/ £2. _l_ °"8 



psca 
-'^oct 



V 



an 



"0 I £8_ 

'VIZ x/6 



(2.105) 



"° / 

V6 / 



and along with singlet state P^ii^g — :^diag((Jo, cxo, (Tq), these define the scalar meson nonet. 
The meson octet matrices are constructed in a similar fashion to the baryon octet matrix; 



Oct 



/ UU U 
U^d 

du 
d ^ s 
\ su 



d ud d — >• s us \ 
dd ds 
sd ss ) 



(2.106) 



The singlet and octet representations of both oj and a (a;o, (^S; '^O; appearing in Eqs. (2.103)- 
(2.105)) are not however the physical particles which we wish to include in the model; these 
are linear combinations of the physical particles. Due to explicit SU(3) flavor-symmetry 
breaking (m^ > m^, m^), a mixture of the unphysical and wq states produces the physical 
w and mesons, while a mixture of the unphysical erg and ctq states produces the physical a 
and /o mesons, the properties of which we list in Appendix B. 

The octet and singlet states are represented by linear combinations of quark-antiquark 
pairs. The state vectors for these are 

l^jg) = ks) = {\uu) + \dd) - 2|ss)) , Iwo) = ko) = ^ (|nn) + \dd)\ + |ss)) , (2.107) 

where the normalizations arise by ensuring that 

(elO = l; {t^8,^0,CT8,CTo}. (2.108) 
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Since the quark contents of the physical states are predominantly 

1 



uj = a 



V2 



{uu + dd), and (j) = fo = —ss, 



(2.109) 



we can replace the octet and singlet combinations with the physical states via the replace- 
ments of 



UJs 



0-8 



1 ^ , -U 2 

73 "^ + 75 



1 ^ _L 2 f 
73 ^ + 76 ^0' 



^0 



1 



(2.110) 
(2.111) 



Using these definitions, we can express the octet and singlet states in terms of the physical 
states 



pvec 
^oct 



/ ^ 4_ _L A. 

P~ 
K*- 



V2 



+ 



_|_ 
18 



K*+ 



\ 



+ ^ 



(2.112) 



V9 J 



psca 
oct 



' V2 VT8 V9 



V 



_|_ o" _|_ h_ 



(2.113) 



Each of the above mesons can interact with a pair of baryons in three possible SU(3) invariant 
ways, which we shall identify as F-style (anti-symmetric), D-style (symmetric) and 5-style 
(singlet). The singlet mesons and (Tq are associated with an S-style coupling, while the 
octet particles are associated with F- and D-style couplings. 

To determine these -F-, L*-, and S-style couplings for the vector and scalar mesons we 
need to calculate the SU(3) invariant Lagrangian density coefficients symbolically for each 
isospin group, with each SU(3) invariant combination^^ given by: 



Tr(BPB) - Tr{BBP) 
Tr(SPoctS) - MBBPoct), 

Tr{BPB) + Tt{BBP) - -Tr{BBP)Ti{P) 

3 

TriBPoctB) + Tr{BBPoct), 

Tr(5S)Tr(P) 

TV(55)Tr(Psing). 



(2.114) 



where we have expanded the meson matrices P according to Eq. (2.104), and we note that 
the octet matrices B, and B are traceless. 



Following the notation of Ref. [32]. 
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Together, these terms can be combined to form an interaction Lagrangian density for 
ah possible SU(3) invariant interactions involving each meson isospin group, with octet and 
singlet coupling coefficients F, D, and S (each defined separately for each isospin group) as 

r^^t = -V2{F [BBP]^ + D[BBP]^} - [BBP]^, (2.115) 

where the remaining numerical factors are introduced for convenience. 

If we evaluate this Lagrangian density by matrix multiplication of B, the octet and singlet 
matrices of Pvcc, and B = B'^j^ in the combinations stated in Eq. (2.115), we can extract 
the coefficients of each baryon- meson vertex in terms of F, D and S factors. These are 
summarized in Table 2.1 for vertices involving the physical vector mesons uj and and (j), for 
pairs of like baryons^^. The summary for the scalar mesons is the same under replacements 
of w —)■ cr, p —)• Sq, and (p — )• /q. 



Table 2.1: F-, D-, and S'-style couplings of like baryon-baryon pairs to vector mesons used in 
these models, according to vertices of type B + P ^ B. The summary for the scalar mesons 
is the same under the replacements of — > o", p^ao, and cp ^ Jq. 
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i {6V2D + 
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-6V2D + V3S) 
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70^0, 



oc -D + F 

i {3D + 9F- V6S) 
i {3V2D + 9V2F + ^35) 
oc D-F 

oc i (3-D + 9F - V6S) 

oc i {3V2D + 9V2F + ^35) 



oc 



oc 



^ As discussed in Section 2.2, any fiavor-changing meson-baryon interactions would produce a null overlap of 
ground-state operators, and as such we only focus on the like-baryon interactions of the form gsa'tpBCitpB' ^bb' 
for a meson a in this discussion. 
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The physical (p and /o states are purely strange quark components. These do not couple 
to nucleons significantly (since nucleons contain only up and down valence quarks) the only 
way to produce these mesons is via gluons. Thus we set the (normalized or not) couplings of 
these mesons to zero; = gsfo — 0- We are then left with the physical a and u mesons as 
the effective meson degrees of freedom. 

If we denote the total (but not normalized) coupling (now including all prefactors of 
Eq. (2.115)) of a (like) baryon pair BB to a meson a by fsa, and we calculate the SU(3) 
invariant combinations for the singlet ujq and mixed state wg, we can use the relation be- 
tween these normalizations from Eq. (2.107) to relate the S'-style couplings to the remaining 
couplings via 

- I = fNuo = V2/JVC.8 = ^{D-3F), (2.116) 

and so we can reduce the relation of the couplings to 

S = V6{3F-D), (2.117) 

and we can therefore find the couplings of the singlet ivq meson in terms of just F and 
D factors. After removing the strange quark components and substituting the result of 
Eq. (2.117) we obtain a summary of couplings as shown in Table 2.2. 

We note however that the couplings in Tables 2.1 and 2.2 do not display isospin symmetry 
manifestly, though our original Lagrangian density (refer to Section 2.1) was constructed in 
terms of isospin groups only with common coefficients. This can be remedied by considering a 
general Lagrangian density constructed from isospin groups, which we shall restrict to terms 
involving like baryons and the mesons we are interested in, to give 

-Cr^t = -fNpiNr^N) ■p + if^pitxt)-p- /Hp(HfTH) • p 

-fNUNN)io - /A^(AA)a; - • S)a; - /=^(HH)^^, (2.118) 

where the A^, A, and H isospin groups are defined as before as 

N=(P], A=(A), H=(r°). (2.119) 



Table 2.2: Couplings of like baryon-baryon pairs to vector mesons used in these models, 
according to vertices of type B + P ^ B using the relation of Eq. (2.117) 
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The p mesons terms are defined in isospin space as linear combinations of the physical charged 
states (as we did in Section 2.3.1) as 

or, equivalently as 

Pi = ^ip+ + p'), P2 = ^ip--P+), P3 = P°, (2.120) 

with the same convention for the replacement of p— > S. This gives the expansion 

S-p = S+p- + (2.121) 
We can expand the Lagrangian density term by term to find the individual interactions 

(Nt'^N).p = {p n)T^p'(^P 

= {pn + np)pi + i(j)n — np)p2 + {pp — nn)p3 

= -^(P^ + '^P)iP+ + P~) - -^{PIT- - np){p~ - p+) + {pp - nn)po 



pppo - nnpo 



+ V2pnp+ + V2npp_ , (2.122) 



where we note that a term BBP indicates the annihilation of a baryon B with a meson P, 
and the creation of a baryon B according to the reaction B + P ^ B. Continuing to expand 
terms, for the S baryons we have 



(SxS)-p = -ip^ (|s°S--S+So) -^pO(|s+E+-S"S- 

(2.123) 



and for the S baryons, 

(HfT-).p = {T^ '^)rjp^(f_ 



+ E~p)p, + i(HOH- - E~E^)p2 + (HOH" - E-E-)p: 



3 



■0-Vo - E~E~po + V^E~p+ + V^E°p_. (2.124) 



The iso-scalar terms are more straightforward; 

(NN)uj = ppoj + nnuj, (2.125) 

AAoj (requires no expansion), (2.126) 

(f-S)a; = S+S+Lj + SOS^w + S^S^cj, (2.127) 

{EE)co = EPE^oj + E^E^uj. (2.128) 
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Once we have calculated the full interaction Lagrangian density, and the F and D coeffi- 
cients of each interaction, we have factors of the following form: 

Ant = y^^y^^ ABmfBmXBrn', XBm = CBmBBm, (2.129) 
B m 

where A^^p = i, and for all other interactions ABm = — 1- The term Xbtu is the expanded 
interaction term (after expanding cross products, etc.) arising from the general Lagrangian 
density, Eq. (2.118), and contains factors of CBm = ±l,ib\/2. We also require a term 
calculated from the SU(3) invariant combinations MBm] the coefficient of the interaction 
BBm in terms of F and D factors as found in Table 2.2. 

To calculate the values of fBm we apply the following formula: 

fBm = ^MBm. (2.130) 
ABm 

For example, consider the interaction vertex w + S*^ — )■ S'^: 

2D , +1, 2D, 2D 

Ae^ = -1, Cs. = +1, A/ec. = -^, hc. = —^{-^) = ^- (2.131) 

Performing these calculations for every possible interactions provides (consistently) the fol- 
lowing couplings of the octet of baryons to the octet of mesons: 



fNp = D + F, /ap = 0, /ep = 2F, hp = F-D, 
fNu = 3F-D, f^^ = -lD + 2F, h^ = 2F, fsu = F - D. 



(2.132) 



We can further simplify our calculations by examining all the different currents that one 
can form using a baryon, an antibaryon and a meson. As discussed in Ref. [33], all possible 
couplings of baryons to vector mesons (denoted by BBV) should be considered when writing 
out the most general Lagrangian density. By calculating the currents (prior to making any 
approximations or assumptions that appear in earlier sections here) we are able to find the 
F-, D-, and S'-style couplings of the form BBX where X is a meson with either scalar (S), 
vector (V), tensor (T), axial- vector (A) or pseudo-scalar (P) spin form. 

Under an expanded SU(6) spin-flavor symmetry, the currents are shown in Table 2.3, 
where the various vector couplings are of the forms 

Vi = y2 = ipcTpuq^i', V3 = i'Qfii', (2.133) 

and we use the convenience definitions of 



i 4M^ -I- a 



If we now consider this as a low energy effective field theory, we can consider the case 
of = 0. We can also enforce rotational symmetry due to lack of a preferred frame (or 
direction) and thus remove the spatial components of both the mesons and the momenta, 
so that = {V^,0) and g'^ = {q^,0), as per Section 2.3.1. Along with cjoo = 0, all terms 
proportional to q'^ vanish, and H = 2. Using these assumptions, the currents are reduced to 
those found in Table 2.4. 
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Table 2.3: F-, D- and 5-style Baryon currents for all types of meson vertices of the form 
BBX where X is a meson with either scalar (S), vector (V), tensor (T), axial (pseudo-) vector 
(A) or pseudo-scalar (P) spin form. Adapted from Ref. [33]. 

F D S 

S \H^i; \Hi}^lj 

V2 -^itpaf^^i; +^iii)a^yi> --^^ii^cFtiui^ 

V3 

A ^HTp-f^jf.il; ^Hijj-f^-ff.il; lH'4}'y5-/^^p 

P lH4)j5ip lHip-/5ilj \Hi>-i^ilj 



We can now observe the relations between the F- and D-style couplings (with the 5-style 
coupling now contributing to F and D)] First, as a check, we observe that the ratio D/F 
for the pseudo-scalars (and the axial- vectors for that matter) is indeed | as commonly noted 
in the literature [34,35] under SU(6) symmetry [36]. Less commonly found in the literature 
is that the 7^-type vector coupling is purely F-style, thus the vector analogy of the above 
relation is D/F = 0, implying D = 0. 

Using the couplings of Table 2.2, we can evaluate the couplings of the vector mesons to 
the entire baryon octet. This provides us with a unified description of the couplings in terms 
of an arbitrary parameter F. These couplings are thus 

fNp = F, fAp = 0, hp = 2F, fsp = F, 

2.135 

fNu. = ^F, /a^ = 2F, h^ = 2F, h^ = F 

We can normalize these results to the nucleon-o; coupling, since we will fit this parameter to 
saturation properties (refer to Section 2.8). Thus the normalized couplings are 

gBm = QNu^r^- (2.136) 

jNui 

We can then separate the meson couplings, since the the normalization above results in the 
following relations, using isospin I^, and strangeness 5"^ of baryon B; 

(3 -5b) 2Ib in ^1-7\ 

QBuj = ^ gNuj, 9Bp = gNuj- [Z.i6() 

These results are consistent with a commonly used naive assumption that the uj meson couples 
to the number of light quarks, and that the p meson couples to isospin. To emphasize the 
isospin symmetry in our models, we will include the isospin as a factor in our Lagrangian 
densities in the form of the r matrices. In doing so, rather than having an independent 
coupling for each isospin group, we will have a global coupling for the p meson, gp. 
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Table 2.4: D- and S'-style Baryon currents with mean-field assumptions = (V^,0), 
qi^ = (gO,0), and = Q. 
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Similarly to the above relations for the vector mesons, we have the same relation for the 
scalar mesons; that the coupling is purely F-style {D = 0). Therefore the couplings for the 
scalar mesons are the same as for the vector mesons, under the replacements a; — )• o", p — t- ao. 
In the calculations that follow, we shall further neglect the contributions from the scalar 
iso- vector oq due to their relatively large mass (refer to Table B.l). 

As an alternative to the SU(6) relations for the p meson coupling gp, we can use an 
experimental constraint. As we have shown above, the p meson couples to isospin, and 
as we will show in Section 3.1 the isospin density is proportional to the asymmetry between 
members of an isospin group; for example the asymmetry between protons and neutrons. This 
asymmetry is measured by the symmetry energy 04 = flsym (derived in Appendix A. 5. 4) which 
appears in the semi-empirical mass formula (the connection is derived in Appendix A. 5. 2) 
which in the absence of charge symmetry is defined by Eq. (2.74). The coupling of p to 
the nucleons is found such that the experimental value of the asymmetry energy of asym = 
32.5 MeV is reproduced at saturation. The coupling of p to the remaining baryons follows 
the relations above. 



3 
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Models Considered 



As outlined in the introduction, QCD is widely believed to be an accurate description of 
strong-interaction particle physics. As a non-perturbative theory, it cannot be solved ana- 
lytically. In order to make any predictions for this theory we must simulate the physics of 
QCD, and we choose to do so using a model. 

As with most fields of research, there are several choices for models to investigate. The 
validity of these models must always be challenged, and there is always a tendency to have a 
preference for a particular model. Our goal is to work with a model which does not introduce 
any physics that is not manifest in Nature, and does not make any assumptions that cannot 
be verified. To this end, we begin with a model called Quantum Hadrodynamics (QHD) 
which — although it is not a quark-level model — simulates the fundamental interactions of 
QCD with effective meson interactions. 

3.1 Quantum Hadrodynamics Model (QHD) 

The origins of the mean-field approximation and QHD reach back to the non-relativistic work 
of Johnson and Teller [37], which was reformulated by Duerr in a relativistic model [-38]. Once 
Chin and Walecka [9,39] successfully reproduced saturation properties (refer to section 2.8), 
QHD as it is known today was born. The formalism for QHD used for this work is expertly 
detailed by Serot and Walecka [18] and by Furnstahl and Serot [40]. QHD was the first great 
step towards a particle-physics understanding of nuclear matter, particularly in the form of 
neutron stars. 

QHD is an effective^, fully relativistic'^ field theory which makes use of a mean-field 
approximation (MFA, refer to Section 2.2) to describe Dirac nucleons interacting at the 
quantum level. The gluons of QCD are simulated by a delicate balance between attractive 
interactions of scalar mesons, and repulsive interactions of vector mesons which, when added 
together, produce an effect which approximates the strong nuclear interaction; the interaction 
responsible for holding protons and neutrons together inside a nucleus where the the Coulomb 
interaction would otherwise cause the protons to repel and prevent any nuclei from existing. 

In the original formulation of QHD (later dubbed QHD-I), interactions between the de- 
generate'^ iso-doublet of nucleons (protons p, and neutrons n) involved the scalar-isoscalar 
a and vector- isoscalar lo mesons in the zero temperature limit. This was soon expanded to 
QHD-II by inclusion of the uncharged vector-isovector po meson. We extend this to include 
the full baryon octet by including the hyperons; A, T,^, S~, and H~. For historical 
purposes, and to emphasise the effects of each of these advances, when we present results we 
shall do so for the most sophisticated version of a particular model, but also for conditions 
corresponding to these variations for comparison. 

We further extend this to use the SU(6) spin-flavor symmetry to relate the couplings of 
all the baryons to all the mesons using the F-style couplings as detailed in Section 2.12. 

^Not involving fundamental particles as the degree of freedom, but rather composite particles which provide 
a useful approximation to the physics. 

^Making use of and obeying relativity, both special and general. 

'^In which the particles share a common mass, and thus satisfy charge symmetry. 
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3. 1. Quantum Hadrodynamics Model (QHD) 



We make our model more physically realistic by using the physical masses for the baryons, 
since charge symmetry is violated in Nature. This provides a mass difference between the 
protons and neutrons, and the model becomes sufficiently precise that we may include leptons, 
which have a considerably smaller mass than the baryons. We may then consider the effects of 
/3-equilibrium between species; by including leptons'* {£ € {e~,/U^}) we allow the possibility 
of considering a balance between various charged species, although we are not explicitly 
modelling charge-conserving interactions, as this would require the inclusion of photon terms 
in the Lagrangian density, photons being the mediators of the electromagnetic force. If 
we did not include the leptons, we would not be able to consider globally charge neutral 
nucleonic matter, as no negative charges would be available to balance the positive charge of 
the protons. 

We will use the term 'configuration of a model' to indicate differences (such as types of 
particles included or neglected) within a particular model. The following discussion will focus 
on the most sophisticated configuration — the octet of baryons in /3-equilibrium with leptons 
— and less sophisticated configurations can be obtained by restricting this description. The 
Lagrangian density that describes such a configuration of QHD is 

k 

t 

where the indices k G {A^, A, S, H} and I G {e~, represent the isospin group of the baryon 
states and the lepton states, respectively, fj-;,) are the isospin matrices for each isospin group 
(refer to Eq. (2.12)), and ^/^^ corresponds to the Dirac spinors for these isospin groups, i.e. 

V'iv = (^;;), ^A = w, V'E=||^soj, V's = (^:°). (3.2) 

The vector field strength tensors are 

F^^ = d^u'^ - d-uj^ = d^Pa - - gpea^cP^P^'c, (3.3) 

V'f is a spinor for the leptons, and 5C are renormalization terms. The values of the baryon 
and meson masses (in vacuum, as used in the calculations herein) are summarized later in 
Table 5.1. We have neglected nonlinear meson terms in this description for comparison pur- 
poses, though it has been shown that the inclusion of non-linear scalar meson terms produces 
a framework consistent with the QMC model without the added hyperfine interaction [41] 
(see Sec. 3.2). 



■'We can safely neglect the contribution of r leptons — which have a mass in excess of 1776 MeV [17] which 
makes them more massive than the H hyperons — since their chemical potential would be equal to that of 
the electrons, and thus one would require an extraordinarily large electron contribution in order to provide a 
non- negligible t Fermi momentum. 
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Assuming that the baryon density is sufficiently large, we use a Mean-Field Approxi- 
mation (MFA, as described in Section 2.2) with physical parameters (breaking charge sym- 
metry) in which the meson fields are replaced by their classical vacuum expectation values, 
a — )• (a)ciassicai- With this condition, the renormalization terms 6C can be neglected. 

By enforcing rotational symmetry (refer to Section 2.3.1) and working in the frame where 
the matter as a whole is at rest, we set all of the three-vector components of the vector meson 
fields to zero, leaving only the temporal components. Furthermore, we remove all charged 
meson states as per the discussion in Section 2.3.1. Consequently, because the mean- fields 
are constant, all meson derivative terms vanish, and thus so do the vector field tensors. The 
only non-zero components of the vector meson mean fields are then the temporal components, 
(w'^) = {uj)6P^ and {p^) = {p)d'^^. Similarly, only the third isospin component of the p meson 
mean-field is non-zero, corresponding to the uncharged po meson. 

The couplings of the mesons to the baryons are found via SU(6) spin- flavor symmetry [32]. 
This produces the following relations for the a and w couplings to each isospin group (and 
hence each baryon B in that isospin group) as per Section 2.12: 

111 111 

2 dNa = 2 9Aa = 2 9T,a = OEa, ^ 9Nuj = ^ 9Auj = ^ ^Ew = ffHw (3.4) 

Using the formalism of Eq. (3.1) with isospin expressed explicitly in the Lagrangian density, 
the couplings of the p meson to the octet baryons are unified, and thus we can calculate the 
coupling of any baryon to either the a, lo or p meson. 

By evaluating the equations of motion from the Euler-Lagrange equations, 



dC „ dC 



we find the mean-field equations for each of the mesons, as well as the baryons. Prior 
to applying the MFA for the mesons, the a equation of motion produces a Klein-Gordon 
equation, while the oj and p equations of motion produce Maxwell equations. This is by 
construction, and these terms can be found in Appendix A. 3. Returning to the use of the 
MFA, the equations of motion for the meson fields are 



B B 

(P) = E4(V^^7%)3V';.> = E^(V'Ir(,)3V'fc) = E4(V'U3B^s)' (3-8) 

k P k P BP 

where the sum over B corresponds to the sum over the octet of baryons, and the sum over 
k corresponds to the sum over isospin groups. I^b is the third component of the isospin of 
baryon B, as found in the diagonal elements of T(fc)3 in Eq. (2.12). (w), (p), and {a) are 
proportional to the conserved baryon density, isospin density and scalar density respectively, 
where the scalar density is calculated self-consistently. 

The Euler-Lagrange equations also provide a Dirac equation for the baryons 

E [^^ - 9Bu.l\uj) - 9pl''hB{p) - Mb + gBa{<T)] i>B = 0, (3.9) 
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3. 1. Quantum Hadrodynamics Model (QHD) 



in which we have inserted the expressions for the self-energies in QHD. The effective mass 
(as defined by the scalar self-energy inserted above, as per Eq. (2.43)) is given by 



MB + J^%=MB-gBa{(j), 



(3.10) 



and the evaluation of this is shown in Fig. 3.1 for the octet of baryons. We note here that 
this definition of the effective mass — being linear in the scalar field — has the possibility of 
being negative for particular values of {a). We will discuss this issue (and a remedy to it) in 
more detail later, but for now we shall acknowledge that this definition of the effective mass 
includes only the first term of many as a linear approximation. 

The baryon chemical potential (Fermi energy) is defined here as the energy associated 
with the Dirac equation Eq. (3.9), which involves the above self-energies as 



^J'B = (-Fb = \I^?B + {^^hY + gBuj{^^) + 9phB{p)- 
The chemical potentials for the leptons are simply 



(3.11) 



(3.12) 



since the leptons do not interact with the mesons. 

The energy density £■ and pressure P for the EOS can be obtained using the relations for 
the energy-momentum tensor for a perfect fluid at rest, where u'* is the four- velocity 

1 



{8 + P) u^u"" + P^^\ 



3 



(3.13) 
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Fig. 3.1: (Color Online) Baryon effective masses M* as defined by the scalar self-energy for 
QHD. Note that it becomes possible that at some densities these effective masses may become 
negative. Negative effective masses in the Dirac equation imply the presence of antibaryons, 
which this model has neglected, so at this point the model becomes unreliable. 
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since = and uqu^ = — 1, where g'^'^ here is the inverse metric tensor having a negative tem- 
poral component, g = diag(— 1, +1, +1, +1), in contrast to the Minkowski metric of Eq. (2.5). 
In accordance with Noether's Theorem, the relation between the energy-momentum tensor 
and the Lagrangian density is 

rA'- = _gM-£ + 9A«^^^. (3.14) 

We find the energy density and pressure for QHD to be the following sum of contributions 
from baryons B, leptons i, and mesons m to be 

s = 

j=B,e,m. 

i=B,e ^ ' a=a,ijj,p 
j=B,l,m 

(||±i) / P£(M2 A + Y. - 



where Jj is the spin of particle i (Jj = ^ V -i G {B,€\) which in this case accounts for 
the availability of both up and down spin-states, and is the Heaviside Step Function (see 
Eq. (2.39)). Note that the pressure arising from the vector mesons is positive, while it is 
negative for the scalar meson. For a full derivation of these terms, refer to Appendix A. 3. 

The expression for the self-consistent scalar field {a) is determined by the derivative of 
the energy density with respect to the effective mass. In the case of QHD this expression is 
given by 

which is solved self-consistently since the effective mass (appearing in the integral) is defined 
in terms of this quantity, as per Eq. (3.10). 

The couplings qn^j and qmw are determined such that symmetric nuclear matter (in which 
Pp = Pn = O.Sptotal) saturates with the appropriate minimum in the (binding) energy per 
baryon, as per Section 2.8. The couplings for QHD which provide a fit to saturated nuclear 
matter are shown in Table 5.2. 

The EOS for QHD can be obtained by finding solutions to Eqs. (3.6)-(3.8) subject to 
charge-neutrality, conservation of a chosen total baryon number, and equivalence of chemical 
potentials. These conditions can be summarized as 

= Ej QiPi 

P = EiBiPi \ ie{p,n,A,S+,EO,S-,HO,H-,e-,^-}. (3.18) 



Once these equations are solved, the energy density and pressure can be calculated. 
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3.2. Quark-Meson Coupling Model (QMC) 



It should be noted that, as with many relativistic models for baryonic matter, once we 
include more than one species of baryon this model eventually predicts the production of 
baryons with negative effective masses at sufficiently high densities (p > 1 fm""^). This is a 
direct result of the linear nature of the effective mass as shown in Eq. (3.10). As the Fermi 
energy (see Eq. (3.11)) approaches zero, the cost associated with producing baryon-antibaryon 
pairs is reduced and at this point the model breaks down. From a more physical point of 
view, as the density rises one would expect that the internal structure of the baryons should 
play a role in the dynamics. Indeed, within the QMC model, the response of the internal 
structure of the baryons to the applied mean-scalar-field ensures that no baryon mass ever 
becomes negative. 

3.2 Quark-Meson Coupling Model (QMC) 

Up to this point we have only considered QHD as an effective field theory; we have considered 
baryons as the fundamental degrees of freedom for this scale. It may however be the case 
that further internal degrees of freedom are more significant at high-densities. We therefore 
wish to extend our model to include the effect of baryon structure in the form of quarks^. 

Deep inelastic scattering experiments have shown that nucleons do indeed have internal 
structure [l.]], and few would discount quarks as the fundamental particles involved. To 
include this degree of freedom, one needs to solve the boundary condition equations for the 
Dirac particles and include the energy contribution from a ' bag' similar to the MIT bag model 
parameter, as will be shown. 

The work here will focus on the latest development of the QMC model [44] which includes 
a quadratic term with numerical factor d in the effective mass which accounts for the scalar 
polarizability. Like QHD, QMC is a relativistic quantum field theory formulated in terms 
of the exchange of scalar and vector mesons. However, in contrast with QHD these mesons 
couple not to structureless baryons but to clusters of confined quarks. As the density of the 
medium grows and the mean-scalar and mean-vector fields grow, the structure of the clusters 
adjusts self-consistently in response to the mean-field coupling. 

While such a model would be extremely complicated to solve in general, it has been 
shown by Guichon et al. [15] that in finite nuclei one should expect the Born-Oppenheimer 
approximation (in which we are able to isolate and distinguish the effects of the quarks, in 
the same way that one is able to separate the net effects of electrons in an atomic calculation) 
to be good at the 3% level. Of course, in nuclear matter it is exact at mean-field level as a 
result of the constant meson fields. 

Within the Born-Oppenheimer approximation, the major effect of including the structure 
of the baryon is that the internal quark wave functions respond in a way that opposes the 
applied scalar field. To a very good approximation this physics is described through the 
'scalar polarizability' d, which in analogy with the electric polarizability^, describes the term 
in the baryon effective mass quadratic in the applied scalar field [47-51]. Recent explicit 
calculations of the equivalent energy functional for the QMC model have demonstrated the 
very natural link between the existence of the scalar polarizability and the many-body forces, 
or equivalently the density dependence, associated with successful, phenomenological forces 

^Some portions of this section are adapted from Carroll et. al. [ I J]. 

"^In QED, the expression for the energy shift due to the quadratic Stark effect is A_E — — ^q|_EzP, in which 
Q is the electric polarizability and Ez is the external electric field, taken to point along the £-axis [46]. 
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of the Skyrme type [52,53]. In nuclear matter, the scalar polarizability is the only effect of 
the internal structure in the mean-field approximation. On the other hand, in finite nuclei 
the variation of the vector field across the hadronic volume also leads to a spin-orbit term in 
the nucleon energy [45]. 

Once one chooses a quark model for the baryons and specifies the quark level meson 
couplings, there are no new parameters associated with introducing any species of baryon 
into the nuclear matter. Given the well known lack of experimental constraints on the forces 
between nucleons and hyperons (let alone hyperons and hyperons) which will be of great 
practical importance as the nuclear density rises above (2-3)pO) this is a particularly attractive 
feature of the QMC approach and it is crucial for our current investigation. Indeed, we point 
to the very exciting recent results [54] of the QMC model — modified to include the effect 
of the scalar field on the hyperfine interaction — which led to A hypernuclei being bound 
in quite good agreement with experiment and S hypernuclei being unbound because of the 
modification of the hyperfine interaction, thus yielding a very natural explanation of this 
observed fact. We note the success that this description has generated for finite nuclei as 
observed in Ref. [53]. 

While we will use the QMC model for our considerations of baryon structure here, we 
note that there has been a parallel development [55] based upon the covariant, chiral symmet- 
ric NJL model [50], with quark confinement modelled using the proper time regularization 
proposed by the Tiibingen group [57,58]. The latter model has many advantages for the com- 
putation of the medium modification of form factors and structure functions, with the results 
for spin structure functions [59, GO] offering a unique opportunity to test the fundamental 
idea of the QMC model experimentally. However, in both models it is the effect of quark 
confinement that leads to a positive polarizability and a natural saturation mechanism. 

Comparisons between the QHD and QMC derivations a posteriori reveal that although the 
underlying physics of QHD and QMC is rather different, at the hadronic level the equations to 
be solved are very similar. Full discussions and derivations can be found in Refs. [41,47,61,62]. 
We shall rather focus on the changes to QHD which are required to produce the QMC model: 

1. Because of the scalar polarizability of the hadrons, which accounts for the self- 
consistent response of the internal quark structure of the baryon to the applied scalar field [53] , 
the effective masses appearing in QMC are non-linear in the mean-scalar field. We write them 
in the general form 



where the weightings u;^, tD^, and the scalar polarizability of the nucleon d, must be cal- 
culated from the underlying quark model. Note now that only the coupling to the nucleons 
gNcT, is required to determine all the effective masses. 

The most recent calculation of these effective masses, including the in-medium dependence 
of the spin dependent hyperfine interaction [51], yields the explicit expressions: 



M*B = Mb - w^B gNa{a) + -w% {gNa{cy)? 



(3.19) 



MN{{cr)) = Mat - gNa{(y) + 0.0022 + 0.1055i?Xr'='= - 0.0178 



(3.20) 




(3.21) 
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Ms((ct)) 



0.6706 - 0.0638/2^^}'='= - 0.008 (r%'"'] 

) 



2 



-0.0007 + 0.0786iZ^^^ - 0.0181 (r%^'') 



9Na{cr) 
{9Na{o)f , (3.22) 



Mh((ct)) 



+ 



0.3395 + 0.02822/?^'='= - 0.0128 (r%'''' 



-0.0014 + 0.0416i2^^^ - 0.0061 (R%''^y 



9Na{cr) 



i9Na{cr)y 



(3.23) 



We take R^^'^ = 0.8 fm as the preferred value of the free nucleon radius, although in practice 
the numerical results depend only very weakly on this parameter [53]. 

Given the parameters in Eqs. (3.20)-(3.23), all the effective masses for the baryon octet 
are entirely determined. They are plotted as functions of the Hartree scalar self-energy 

= —gNa{<7) in Fig. 3.2 and we see clearly that they never become negative (note that the 
range of covered here corresponds to densities well above (6-8)/9o in QMC). 




Fig. 3.2: (Color Online) Baryon effective masses for QMC as a function of the QHD (linear) 
scalar self-energies. The values at S?^^ = are the vacuum masses as found in Table 5.1. In 
order to emphasise the curve we have shown the effective masses beyond = 2000 MeV, 

though in practice we only require values up to — S?^^ = 800 MeV which corresponds to 
densities of ~ 2 fm~'^ (6-8 po)) beyond which higher order terms not shown in Eq. (3.19) 
become significant. The (shaded) unphysical region of S?^^ > contains a point for which 
all of the effective masses become unified; this is the point at which SU(3) symmetry is an 
accurate symmetry. 
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2. The mean-scalar field {a) is derived self-consistently by taking the derivative of the 
energy density with respect to (cr), thus the scalar field equation 

= E ^C((.))ef±ll / Mimi^ A, (3.24) 



has an extra factor when compared to Eq. (3.17), denoted by 

C{{a)) = wl- w^dgNaicy). (3.25) 

Note that the the scalar polarizabihty d term in C{{a)) does not have the factor of ^ that is 
found Eq. (3.19), because of the differentiation of the squared term. 

Given this new term in the equation for the mean-scalar field, we can see that this allows 
feedback of the scalar field which is modelling the internal degrees of freedom of the baryons. 
This feedback prevents certain values of {a) from being accessed. 

3. The couplings to the nucleons are re-determined by the fit to saturation properties (as 
per Section 2.8) with the new effective masses for the proton and neutron. The couplings for 
QMC which provide a fit to saturated nuclear matter are shown in Table 5.2. 

Given these changes alone, QHD is transformed into QMC. The implications of these 
changes however will be profound. 



3.3 MIT Bag Model 

We consider two models for a deconfined quark matter phase, both of which model free 
quarks. The first model, the MIT bag model [(i-'5], is commonly used to describe quark matter 
(and, by extension, hadronic matter) because of its simplicity. 

In this model we construct baryons as three quarks confined within a 'bag,' a region defined 
by the local energy density which is greater than the energy density of the surrounding vacuum 
by a factor of B — which is commonly quoted to a power of one-quarter, thus B'^/^ ~ 180 MeV 
— and which is lower in pressure by the same factor. The motivation for this is that the quarks 
reside in a high-energy, low pressure region which confines them. 

This is the observed property of confinement — that no quarks have ever been observed, 
either directly or indirectly, in isolation — that perturbative approaches to QCD are unable 
to reproduce. This property is only of concern for hadronic matter though; at high energies, 
QCD predicts that the quarks become 'asymptotically free', and at this point we can consider 
quark matter. 

The simplest Lagrangian density that can be constructed for such a scenario is 

>Cmit = i'q [i^ - mq) Tpq + B, (3.26) 

where ipq is a spinor for quarks with mass ruq where we typically use masses of ruu = 3 MeV, 
rud = 7 MeV, and = 95 MeV for the up, down, and strange quarks respectively, and B is 
the aforementioned bag energy density. In a similar fashion to the derivations for QHD and 
QMC, the energy density and pressure can be calculated to be 

^ = ^ + E ^^7^^ / ^(^^. - I^dV^^+K? d'k , (3.27) 
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P^.B^y:!2MJ^[llMzBA. (3.28) 

In order to utilise statistical mechanics and explore phase transitions involving this sce- 
nario, we model point objects with fixed masses which possess chemical potentials related to 
the independent chemical potentials of Eq. (2.45) via 

12 11 , , 

f^u = g/in - 2^e, fJ'd = gA^n + ^/^e, fJ-s = fJ-d, (3.29) 

where quarks have a baryon charge of ^ since baryons contain 3 quarks. Because the current 
quarks do not interact with mesons in this model, the quark chemical potential has the same 
form as the lepton chemical potential (no meson terms) and thus 



M<? = Y ^F, + "^q ; q€{u,d,s}. (3.30) 

The EOS can therefore be solved under the conditions of Eq. (3.18). In our calculations 
of this model, we use the current quark masses as found in Ref. [17], which represent the 
physical, bare quark masses. 



3.4 Nambu-Jona-Lasinio Model (NJL) 

As an alternative model for deconfined quark matter, we consider the Nambu-Jona-Lasinio 
(NJL) model [56], in which the quarks have dynamically generated masses, ranging from 
constituent quark masses at low densities to current quark masses at high densities. 

Dynamical breaking of the chiral symmetry produces Nambu-Goldstone bosons — the 
triplet of pions in the two-flavor case — which are massless while the symmetry is preserved, 
but which possess a non-zero mass when the symmetry is broken. According to the Gell- 
Mann-Oaks-Renner relation (for example, see Ref. [64]), in the leading order of the chiral 
expansion, 

"^1 = ^2 {'^u{uu) + md{dd)) , (3.31) 

In 

or more commonly reduced to the fact that the quark mass scales as the square of the pion 
mass, 

niq ~ m^. (3.32) 

At large densities, manifest chiral symmetry is expected to be partially restored. If it was 
fully restored, the quarks would be massless (in which case, jn^r = iriq = 0). The symmetry is 
however not exact, and the quarks retain a very small mass. This small mass is the current 
quark mass as noted in Ref. [17], whereas the dynamically generated quark masses under 
broken chiral symmetry are the constituent masses that together naively sum to the mass 
of a baryon. By using the NJL model rather than a simple treatment for the quark masses 
we endeavor to make our phase transition models more realistic and more sophisticated by 
including more physics believed to represent Nature. The NJL model is a simple construction 
that displays the correct phenomenology, namely DCSB. 
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If we consider a particular choice for a massless Lagrangian density' with strong couphng 
G, in a Lorentz-covariant frame (with no matter or background fields) to be 

C = -^i^^l, + G [{i^i^f - (^75^)'] (3.33) 

we can see that this is invariant under the vector and axial vector symmetries of Eqs. (2.59)- 
(2.60). For example, 

= (3.34) 

where we recall Eq. (2.61), and that 'ip = V'^To- 

The self-energy of this particular model is shown in Fig. 3.3 in which the four-fermion 
vertex is the 1-PI vertex (refer to Section 2.7), and the loop is the full quark propagator. 

If we then approximate the four-fermion vertex to be the bare vertex with coupling G, 
and the loop to correspond to the bare propagator, then the loop now corresponds to the 
quark condensate {ipip) and the self-energy is given by the Feynman diagram as shown in 
Fig. 3.4. The Lagrangian density for NJL then becomes 

C = -i) {i^ + 2G{ti;tl;)) ^ = -i) {i^ + S^) i). (3.35) 

where we can see from Fig. 3.4 that the self-energy is a Dirac-scalar term. We can now 
identify the dynamical quark mass as a scalar self-energy term in this Lagrangian density. 

If we include a further constant mass term tpruQij: in Eq. (3.33) that explicitly breaks the 
chiral symmetry, we can define the effective (dynamic) quark mass as 

m* = ?no -|- Ss = mo — 2G(^^/^). (3.36) 

This is called the 'gap equation' in analogy to superconductivity*^. 




Fig. 3.3: Feynman diagram for the NJL quark self-energy, where the four-fermion vertex is 
1-PI, and the loop corresponds to the full quark propagator. 

^Following the considerations of Nambu and Jona-Lasinio [56]. 

*The BCS theory of superconductivity [65] is a primary motivation for this model, in which electrons in 
a metal may become a paired bosonic state — a Cooper pair — by possessing a lower energy than the Fermi 
energy. In that case, a temperature-dependent energy gap exists, and electron excitations must be of a 
minimum energy, as opposed to the continuous spectrum that the electrons would normally have. 
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Fig. 3.4: Feynman diagram for the NJL quark self-energy, where the four-fermion vertex is 
the bare vertex with couphng G, and the loop corresponds to quark condensate. 



To connect with Section 2.7, we can understand the NJL self-energy of Fig. 3.4 in analogy 
to that of Fig. 2.2 in the (albeit, unphysical) case that we shrink the gluon loop in Fig. 2.2 to a 
point, and in which case the Lagrangian densities are related via the (arbitrary) replacement 
of 

{tpTijf 2i;Ti; {t/jTt/j), (3.37) 

where F is any gamma matrix appearing in the interaction Lagrangian density. The term 
{'4)TiIj) denotes the ground-state (vacuum) expectation value of '4jTiI^, and as noted in Sec- 
tion 2.3.2, the ground-state is a parity eigenstate, thus the only term in Eq. (3.33) that has 
a non-zero vacuum expectation value is {ipip). 

The expression for the Feynman diagram in Fig. 3.4 is simply 

S(p) = 2G{ijij) = 2iGTrSF{0) (3.38) 

where SpiO) = Sf{x — y)5{x — y) is the Feynman propagator for a quark (refer to Ap- 
pendix A. 2) which starts and ends at the same space-time point. 

As with Section 2.7, we have described the above model for a Lorentz-covariant frame, 
in which there are no matter or background fields. In this case, the Feynman propagator in 
Eq. (3.38) is the free-space propagator. We make a further approximation by replacing this 
propagator with the in-medium propagator, one which is in the presence of matter fields, and 
thus we introduce a Fermi momentum. The expression for the quark condensate in this case 
b GconiGS 

K fmie{kF-\k\)9{A-kF) 



ii^gijg) = iTrSFiO) = -4 — ^ / " /^^^ (^-^9) 

^ ^ ) J W A-2 _|_ (j7T,*)2 



where we have introduced a momentum cutoff of kp < A to regularize this integral, at which 
point we expect to recover current quark masses, and Mc is the number of color degrees 
of freedom of quarks. In order to calculate the effective quark mass at each density, we 
must first find the coupling G which yields the appropriate constituent quark masses in free 
space {kp = 0). The coupling is assumed to remain constant as the density rises. Given 
a (free-space) constituent quark mass m^'^^, we can solve Eqs. (3.36) and (3.39) to find the 
coupling 



^ ^ (mf''"'^ - mp) 



Mc fml'^^e{\k\-kF)e{A-kF). 



(27r)3 



d^k 



P + (mfrcc)2 



(3.40) 
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We evaluate Eq. (3.40) for Afc = 3 to produce constituent quark masses of m^'^ = 350 MeV 

using current quark masses of m^''^ = 10 MeV for the light quarks, and a constituent quark 
mass of m^J^^ = 450 MeV using a current quark mass of = 160 MeV for the strange quark 
(both with a momentum cutoff of A = 1 GeV) as per the phenomenology of this field. At 
kp = we find the couplings to be 

Gu,d = 0.148 fm^, Gs = 0.105 fm^. (3.41) 

We can now use these parameters to evaluate the dynamic quark mass m* for varying val- 
ues of kp by solving Eqs. (3.36) and (3.39) self-consistently. The resulting density dependence 
of m* is illustrated in Fig. 3.5 where we observe that the quark masses — particularly the light 
quark masses — eventually saturate and are somewhat constant above a certain Fermi mo- 
mentum (hence, density). For densities corresponding to kp > A, the dynamic quark mass 
is constant. 

We can now construct the quark matter EOS in the same way as we did for the MIT bag 
model, but with density-dependent masses rather than fixed masses. 




J 1 I I I I I 



0.0 0.2 0.4 0.6 0.8 1.0 

K [GeV] 

Fig. 3.5: (Color Online) Dynamic quark masses in the NJL model. The mass at kp = is the 
constituent quark mass, and the mass at the cutoff oi kp = A = 1 GeV is the current quark 
mass. This model successfully reproduces the behavior found within the Schwinger-Dyson 
formalism for dynamical chiral symmetry breaking. 
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3.5 Fock Terms 

As an extension to the models described above, we now introduce terms which are of a higher 
order — the Fock terms — leading us to a Hartree-Fock description of matter. As we wish to 
work to the next leading order, we require a perturbative description of the models in which 
we can identify the next terms. To this end, we utilise Feynman diagrams to illustrate which 
terms shall be calculated, and which shall be neglected, since we can identify the order of a 
diagram by the number of baryon-meson vertices, each of which contributes a factor of the 
baryon-meson coupling. 

The derivation of QHD as described in Appendix A. 3 does not lend itself to such an 
expansion, since all the contributing terms must appear in the Lagrangian density ab initio, 
but we can re-derive QHD perturbatively which will allow the identification of the next-order 
terms. The diagrammatic derivation for QHD at Hartree level can be found in Appendix A. 6. 

The first Fock terms that we wish to calculate are the next-to-leading order contributions 
to the baryon self-energies Ss(fc), which as we shall see are momentum dependent. The full 
derivation of this Fock contribution can be found in Appendix A. 7. In order to calculate 
the Hartree-Fock EOS, we consider Dyson's Equation which self-consistently relates the 
full momentum-dependent baryon Green's function (propagator) G{k) to the bare (vacuum) 
baryon propagator G^{k) and the self-energy, as 

G{k) = G^{k) + G^{k)J:{k)G{k). (3.42) 

This can be represented with Feynman diagrams as shown in Fig. 3.6. 

We can write the as-yet undefined self-energy as a sum of terms by defining components 
of the self-energy proportional to the identity matrix I, or gamma matrices and 7 for scalar 
(S*), temporal (S''), and vector (S* = T,^) self-energies, such that the full self-energy (in the 
case where the following terms are defined in-medium) is expanded as 

E(fe) = S^(fc) - 7^S'^(A:) (3.43) 
= S^(|^|,A:0)-7oS°(|fe|,feO) + 7-^S''(|fcU-0). (3.44) 

Note the differences between Eq. (3.43) and Eq. (2.67); the components of the self-energy in 
the Lorentz-covariant case (free-space). In this in-medium case, we now have dependence on 
k, which will lead to a dependence on Fermi momentum. 

Dyson's Equation can be solved formally to give the baryon propagator 



[G{k)r' = -f^ikf' + S'^(A:)) - [M + S^(A:)]. (3.45) 

For convenience, we define the following quantities; 

M*{k) = M + T.\k), (3.46) 

k* = k + T.\k), (3.47) 

E*{k) = ^{k*Y + (M*)2, (3.48) 

r'^ = k^" + T.^'{k) = [k'^ + YP{k),r]. (3.49) 
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G{k) 




Fig. 3.6: Feynman diagram for Dyson's Equation, Eq. (3.42), providing a self-consistent 
description for the full baryon propagator (double line) involving the bare propagator (single 
line) and self-energy (filled blob) . 



We can express the solution to Dyson's Equation as a sum of Dirac and Fermi components, 
where the Dirac contribution is responsible for Pauli blocking, and the Fermi contribution 
accounts for antibaryons [In], which shall be neglected here'' 



G{k) 
Gpik) 
Goik) 



GFik)+GD{ky, 

[Y'kl + M*{k)] {k*^'kl - {M*{k)f + ie) 



[rkl + M*{k)] 



ITT 



E*{k) 



5{kQ-E{k))e{kF-\k\), 



(3.50) 
(3.51) 
(3.52) 



where the energy E{k) is the self-consistent single-particle energy, calculated 'on-shell'^'^ 

E{k) = [E*{k)-TP{k)]k,=E{ky (3.53) 

In the Hartree case, the tadpole diagrams define the self-energy. If we now include the 
exchange contributions, the diagrams for which (for the vector meson terms) are shown in 
Fig. 3.7, the expressions for the total (including all Lorentz forms) a and ui contributions to 
the self-energy become 



S(fc) 

^Ba{k) 



^Ba{k) + T.Buj{k), 



{k) 



IQBu 



igBu 



(2^ 



Tr[G(g)] 



+ 



9Ba G{q) 



B' 



- 1)1. -ml + ie 



Tr[7^G(<7)] , 9Bul^.g^''G{q)^, 



+ 



B' 



mf. 



{k - q)] 



mi + te 



(3.54) 
(3.55) 
(3.56) 



where for each case in the sum, G refers to the baryon propagator of baryon B. Here g^'^ is 
the Minkowski metric tensor (a.k.a. r/"^ in Eq. (2.5)). 



In their absence we will use a momentum cutoff to regularize the integrals. 
^In which case p^p** = M^, and thus p" = E{p). 
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Fig. 3.7: Summation of Feynman diagrams for all possible (vector) interaction terms con- 
tributing to the self-energy. The first row of diagrams are the tadpole terms, and the second 
row are the exchange terms, where in both cases, each term is derived by inserting the full 
baryon Green's function into the previous term, which self-consistently includes the self- 
energy. Similar diagrams exist for the scalar meson interactions. 

For simplicity, we will neglect the contribution from the p meson in these Fock calculations. 
The self-energy can now be separated into terms proportional to those in Eq. (3.43). 



Since we will be using a momentum cutoff to regularize the integrals, we^^ drop the 
antibaryon components of the baryon propagator GF{k), leaving Goik) as the full propagator. 
Some of the integrals can then be performed to produce the mean-field results plus the Fock 
additions 



^%ik,E{k)) 



^Uk,E{k)) 



B' 

1 

4'K'^k 



' a q- 



ml JO 



^Uk,Eik)) = QB. 



B' 



, M*Jq) 
qdq—^i^- 
E*^{q) 

-{2JB' + l)gB 



\9BaQs{k, q) - g%^@Uk, q) 



An'^k 



1 



47r2A;2 



(27r)2 



q dq 



q dq 



d^q 



\gBa^<T{k,q) + ^gluj&Uk,q) 



i;9Bu'^Ak,q)+gBu>'^Lo{k, 



Following the procedure of Ref. [18] 



(3.57) 



(3.58) 
(3.59) 
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with the definitions (for convenience) of 



In 



Ai{k,q)+2kq 



2kq 



1 



Mf^,q) = ■^Ai{k,q)Q^{k, 



Ai{k,q) 



k' + q'+mt-[E{q)-E{k)r, 



(3.60) 

(3.61) 
(3.62) 



for which q = \q\, k = \k\. 

All of these self-energies are evaluated on-shell at the self-consistent single-particle energies 
defined in Eq. (3.53). To further simplify our calculations, we will use the approximation 

= 0, since the momentum-dependence of this term is a power weaker than and S^. 
The only consequence of this is that in Eq. (3.47), k* — )■ k. 

These contributions to the self-energy will affect the effective masses of the baryons in 
QHD as additional terms in Eq. (3.10), and once we have made this change we can once 
again calculate the energy density. In this perturbative description of QHD, although the 
energy density is still derived using the energy- momentum tensor, the expression for this now 
involves the meson propagators. Using the Hartree meson propagators (since these are not 
affected by the Fock addition to the self-energy) we can evaluate the energy density, and we 
must once again re-fit the baryon-meson couplings such that the saturation properties are 
reproduced. The couplings which provide this are given in Table 5.2 for this case of Fock 
term additions to the self-energy. 

The second Fock terms we can calculate are additional terms in the meson propagators — 
the medium polarisation — which will alter the definition of the energy density. The medium 
polarization is named in analogy to the vacuum polarization, the simplest example of which 
is in Quantum Electrodynamics (QED); in vacuum, a photon can spontaneously create a 
virtual^^ particle-antiparticle pair — namely an electron and a positron — which annihilate 
back to a photon, as shown in Fig. 3.8. While they temporarily exist, the charged pair acts 
as an electric dipole and can partially screen an external electromagnetic field. Similarly, the 
medium polarization is a contribution from the creation of a baryon loop which affects the 
matter field. 




Fig. 3.8: Feynman diagram for vacuum polarization in QED. Here the photon creates two 
fermions which are off-shell (virtual particles) which annihilate back to a photon. The charged 
pair acts as an electric dipole, and thus are able to partially screen an external electromagnetic 
field. 



These particles are 'off-shell', meaning that p^p'^ 7^ m^. 
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The three terms contributing to the scalar propagator are shown in Fig. 3.9. The first of 
these is the bare meson propagator. The second is the Hartree contribution, essentially the 
square of a tadpole term. The last of these is the medium polarization, and the scalar and 
vector contributions to this term are given by 

A'{k) = A°n^(A;)A°(A;) (3.63) 

U^{k) = -lY^gl^ j ^^TT[GD{q)GD{k + q)] (3.64) 



ni"(^^) 



Dl^Ii^J{k)DUk) 

_dS_ 
(2vr)4 



B 



9bu 



Tr 



(3.65) 
(3.66) 



where Ho- and 11;^'^ are the medium polarizations. 

By evaluating the temporal components of the energy-momentum tensor T^^ with the 
Hartree-Fock propagator as shown in Fig. 3.9 we find the full expression for the Hartree- 
Fock energy density 



E 



(2Jb + 1) [^^ 



B 



+ 



1 {2Jb + 1) 



(fkcP'q 



2 (2^)6 70 70 E*{k)E*{q) 
X |5L^0(fe-g) ]^-[E{k)-E{q)fDl{k-q) [k^^q^ + M* {k)M* {q)] 



[E{k)-E{q)fDl{k 



[k*^'q* - 2M*{k)M*{q)] 



(3.67) 



The first line of terms in Eq. (3.67) are of the same form as the Hartree contributions, but 
these terms are now defined using the second-order self-energies. The double integral arises 
from the medium polarization terms. Eq. (3.67) is fully derived in Section A. 7. 

For consistency, once we include Fock terms to the energy density, we must once again 
recalculate the couplings to fit saturation properties. The values that best reproduce the 
data for this Hartree-Fock description of QHD are summarized in Table 5.2. 

While we have presented the results for the a and uj mesons, it should be noted that other 
meson contributions are possible. We could add the contribution from the p meson with the 
appropriate coupling factors, and furthermore since the Fock terms are of a higher order than 
mean-field, we could now include the pion contribution, and all pseudo-scalars not otherwise 
excluded. For now we shall continue to neglect these terms for simplicity. 

An alternative additional contribution to the energy density is the second-order exchange 
contribution as per Ref. [66] which shall be neglected here for simplicity. 
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+ 




Fig. 3.9: The three leading terms contributing to the full meson propagator (double, dashed 
line) in this case of the scalar meson; the bare propagator A*^ (single, dashed line), the 
(squared) tadpole diagram, and the medium polarization diagram. Similar terms can be 
written for the vector meson also. 
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Methods of Calculation 



As is often the case in theoretical physics, for many of the quantities mentioned in the 
preceding chapters we lack all the information required to evaluate them. Although we may 
know the vacuum mass constants and couplings, there are several quantities which are self- 
consistently defined — where the quantity of interest appears on both the left- and right-hand 
sides of an equation — for example the equation for the effective baryon mass in QHD, 

Mb = Mb - gBa{cr) = Mb - qbu V — o" / — ?TT^ ; — d-'k , (4.1) 

^ J (2^)' ^P + Mf, 

is defined self-consistently^ . Technically, Eq. (4.1) describes Nb coupled, self-consistent equa- 
tions, where Nb is the number of baryons in the model considered. Furthermore, although 
we use the total baryon density ptotai as the control parameter, for non-trivial cases the in- 
tegration limit of kpg is reliant on knowing the density of each species separately, which in 
turn can be reliant on via the chemical potentials (refer to Eq. (2.44)). Equations of this 
form cannot be solved analytically, so we turn our attention to numerical calculations. 

All of the calculations performed in this thesis have been computed in Fortran 90 using 
methods inspired by 'Numerical Recipes in Fortran 90' [67]. Rather than publish the partic- 
ular code used to perform the exact calculations of this work, in the following sections we will 
briefly outline some of the challenges that have been met in performing these calculations, 
and the major steps required to reproduce the simulations. 



4.1 Newton's Method 

One of the most widely renowned methods for iteratively solving non-linear equations is the 
Newton-Raphson method, also known simply as 'Newton's Method'. It has the additional 
benefit that it is suitable for solving self-consistent equations, as above. 

Newton's method is elegant in its simplicity; In order to find the value x* which is a 
solution of f{x) = where / is a general function, we begin by considering the gradient /' of 
a tangent to / at a point xi that is sufficiently close to x* (within the radius of convergence). 
The radius of convergence is defined in this case by the region in which the tangent to / 
intersects the x-axis at a point Xj+i such that \x* — < |x* — as shown in Fig. 4.1. 

The gradient of the tangent is given by 

=^^= (4.2) 

Ax Xi - Xj+i 

where we have used the fact that the tangent to / intercepts the x-axis at Xj+i, and thus 
/(xj+i) = 0. Rearranging this gives the iterative procedure for Newton's method; 

x.+i = x.--^. (4.3) 



^In that the effective masses of all the baryons are required in order to calculate the effective mass of a 
single baryon. 
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Fig. 4.1: (Color Online) Diagram for Newton's method as described in the text. Here, the 
starting value Xi lies within the radius of convergence which for this example is bounded 
somewhere between the two extrema. For example, a tangent at or beyond these points will 
have an x-intercept of Xj+i that is further from x* than Xi is. 



If the function is convergent, for some value of i it must be true that — Xi\ < e for some 
tolerance e, at which point the solution x* is considered to be found. To make use of this 
method, we need to define our functions, and the equations that we wish to find the solutions 
to. 

For QHD, QMC, and in general, we use the total baryon density as the control parameter, 
allowing us to calculate the properties of matter at a specific density. At a given density, the 
total baryon density is therefore conserved (constant), and equal to the sum of individual 
baryon densities; 

Ptotai = y^Pg- (4-4) 

B 

We also require that the total charge of the system, defined by 

Qtotal = ^ QbPb (4.5) 
B 

is constant, and vanishing. With two constraining equations we can use Newton's method to 
find the values of (at most) two quantities constrained by these equations. For the purpose 
of simplicity, we use the proton and neutron Fermi momenta kp^ and kp^ since these appear 
in the integration limits of many quantities. 

In keeping with the discussion of chemical potentials in Section 2.5, we could equivalently 
use the neutron and electron chemical potentials; the choice is arbitrary, though using only 
the nucleon Fermi momenta allows us to neglect conservation of total charge, such as the 
case of nuclear matter where we only need a single constraining equation'^, without loss of 
generality. 



■^The nucleon Fermi momenta are constrained such that the densities are equal in that case. 
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Newton's method is a root-solving method, so our constraint equations need to be written 
in a form of f{x) = 0. This is achieved simply as 

/i(^F„) = Ptotai - X]^^ = ^' (^-6) 

B 

f2ikF,) = Qtot.i -^QbPb = 0. (4.7) 

B 

Lastly, in order to solve the self-consistent equations for the baryon effective masses in QHD 
and QMC, we require an additional quantity that we can allow to vary. Considering QHD, 
since (at Hartree level at least) the scalar field {a) is independent of baryon species (it is a 
sum of contributions from all baryons) we use this quantity in a further constraint equation 

= (.) - E ^ / P;';^)^^;''^^.-!^"!) A = 0. (4.8, 

Note that since pb = kp^/Sir'^, each of the above constraint equations contains kp^ and kp^^, 
so these equations are not only non-linear, but are all also highly correlated. 

Although the above derivation of Newton's method is for a function of a single variable, 
we can extend this concept to an arbitrary number of m functions in n variables by replacing 
the derivative by a Jacobian 





dfil 


dxi 


dXn 


dfm 


dfm 


. dxi 


dXn- 



In order to include the hyperon contributions to these calculations, we do not require any 
further constraint equations, since all of the hyperon Fermi momenta can be derived from that 
of the nucleons, by rearranging Eq. (2.46). The constraint equations above, Eqs. (4.6)-(4.8) 
do however change as they now include terms for the hyperons in the sums. 

We can use the same method to perform the quark matter calculations in this thesis, 
with only slight changes to the constraint equations. The conserved baryon density is now 
defined for quark matter as per Eq. (2.93), and the total charge remains zero, though the 
individual charges Qi now reflect the quark charges; Qu = +2/3, Qd = Qs = —1/3. For 
the calculations of the MIT Bag Model, only these two constraint equations are required, 
since the quark masses are constant. For the NJL Model calculations, the quark masses are 
calculated self-consistently in the same manner as the baryon effective masses above, so the 
equation for a quark condensate, Eq. (3.39) can take the place of the scalar field in Eq. (4.8). 

This method is used to produce all of the calculations for this work at Hartree level, but 
it is not suitable for calculations involving the Fock terms; the additional factors that are 
included at Hartree-Fock level alter the self-consistent equations for the effective masses (via 
the self-energy) such that the effective masses become baryon- and momentum-dependent. 
Because of this, we can no longer rely on the mean-scalar field to produce a constraint 
equation, since we would require 2 + nfc xub equations (n^ = number of momentum samples, 
ub = number of baryons in the calculation). Clearly this is not a feasible route to take, since 
for a reasonable number of momentum samples, say nt = 100 we would require 202 equations 
to just calculate properties of nucleonic matter. 
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The alternative is to find a new method of evaluating our self-consistent equations. The 
method we have used is described next. 



4.2 Steffensen's Method 



Until now we have described solving the self-consistent equations within Newton's method, 
but as stated above the Fock terms would introduce too many new equations for this to be 
workable. To make matters worse, the addition of the Fock terms makes the effective masses 
momentum dependent, and so the equations for these become integral equations 



M*B{k) = Mb -93 



gsa f^B' (2JB, + l)M*,(g) ^3 



E 



B' 



m^a Jo 



(27r)3 E*^,{q) 



^9L0<t(^, q) - 9Bw®Uk, q) 



dq. 



(4.10) 



where the terms in this equation are defined in Section 3.5. 

As a solution to this problem, we find a new method of solving the self-consistency. Many 
improvements to Newton's method are available, and we have chosen to investigate one in 
particular; Steffensen's Method. The advantage that this method has over Newton's method 
is that it doesn't require a calculation of the derivative of the function /. 

As an improvement over Newton's method, we replace the derivative (the slope of the 
tangent to /) in the denominator of Eq. (4.3) with the slope g{x) of a line joining the point 
(xj, f{xi)) with an auxiliary point (xj -|- Jxj, /(xj -|- 5xi)) as shown in Fig. 4.2. The choice of 
6xi is somewhat arbitrary, and the choice of 5xi = f{xi) merely ensures that the step size is 
scaled appropriately according to the distance from the solution, and in fact 6xi vanishes as 
f{xi) approaches 0. 




Fig. 4.2: (Color Online) Diagram for Steffensen's method as described in the text. This 
method has the numerical advantage that it does not require the evaluation of a derivative, 
though the requirement of a 'good' starting value is now somewhat more strict. 
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The slope of this new hne is then 

g{x) = /(^^ + /(^^))-/(^^). (4.11) 

f[Xi) 

If we insert this into Eq. (4.3) in place of the derivative we obtain the iterative procedure 
for Steffensen's Method, where the x-intercept of g{x) defines the next approximation to 
fix*) = 0; 

fjXi) fjXjf 

9{xi) f{xi + f[Xi)) - f{xi) 

This method can be applied to solving integral equations such as Eq. (4.10). 

An issue that the reader may notice is that the effective masses in that equation have 
a momentum-dependence in the form of a 1/k term, as well as fc-dependence in the Q{k,q) 
terms; this makes the definition of M* (0) very critical, since the /c-dependence of the term 
(as defined in Section 3.5) produces Qi{0,q) = ln(l) = 0, thus we have a 0/0 issue. Since 
we require the value of M*(0) in order to calculate the integrals above, we must define this 
value very carefully. Since the function is fairly smooth, we use a linear interpolation of the 
function evaluated at the next two momentum points to define the limit of the function at 
k = 0. 

This gives us all the information we need in order to solve the self-consistent equations, 
but we can improve the method by accelerating the convergence, as we shall show in the next 
section. 



4.2.1 Aitken's Process 

We can further accelerate the rate of convergence for Steffensen's Method by using Aitken's 
delta-squared process. If the function we are attempting to find the roots of is f{x) we can 
iterate Steffensen's method to produce intermediate values ho, hi, and h2 such that 



ho = Xi, 
hi = ho 



fjho) 
f{hi 



h2 = hi-^j^. (4.13) 

With these definitions, Aitken's delta-squared process defines the next approximation to the 
root of the function / as 

{h2ho - hif 

Xi+i = ho- . (4.14) 

/i2 - 2hi + ho 

If the process has failed to converge at this stage by satisfying the condition \xij^i — Xi\ < e 
then the entire procedure is repeated for the next value oi x. If convergence is achieved, it is 
trivial to confirm that the convergence has produced the correct value of /(xj+i) = 0. 

Once again, this description has involved only a function of a single variable. The method 
can be modified to handle functions of several variables, but we shall not discuss this modifi- 
cation here as we will only be using this method for nucleon calculations. With this method 
we are able to solve (at least) two self-consistent integral equations, and we shall do so for 
the nucleon effective masses for a suitable number of momentum points. 
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We nest the solving of the self-consistent effective masses within the Newton's method in 
the previous section which still varies the proton and neutron Fermi momenta to find the roots 
of the equations. The functions / that we wish to find the roots of via Steffensen's Method are 
simply rearrangements of Eq. (4.10) for each of the baryons such that /(M*, M*) = 0. Once 
we have solved these self-consistencies, we no longer require the third constraint equation 
Eq. (4.8) in Newton's method. 

4.2.2 Example Integral Equations 

As an example, we consider the fixed-point integral equation 

u{x) = F{u{x)) = + e"^ / u{t) dt, (4.15) 

Jo 

which has the solution 

u*{x) = e'' + l. (4.16) 

Being such a simple case, we can solve Eq. (4.15) in several ways. For the sake of comparison, 
we will compare a naive method; in which a previous function value is used as the next 
guess, the iterative procedure for which is Ui-^i{x) = F(ui), where for this example the 
minimized function is the difference between the left- and right-hand sides of Eq. (4.15), viz 
f{xi) = u{xi) — F{u{xi)); with the Aitken-improved Steffensen's method used for our actual 
calculations. 

In both cases, we supply an identical initial guess of uo{x) = 0, sampled over 1000 x- 
values, and set a tolerance for convergence (via the relative error) at each sampled value of 
Xj (where j = 1, . . . , 1000) of 

^ Ui+l{xj)-U,{xj) ^ ^ ^ 
Ui+l{Xj) 

In both cases the procedures converge to within the above tolerance, and the actual errors 
relative to u*{x) are shown in Fig. 4.3. The naive method converges with 20 iterations, while 
Steffensen's method converges with just 11, a vast improvement. 

The speed of convergence is one advantage of the Aitken-improved Steffensen's method, 
but the more important advantage comes when we try to solve a more complicated integral 
equation, such as Eq. (4.10) for a single baryon, which is a fixed-point integral equation. 

With all constants defined the integral equation reduces to M* = f{M*,kp), where the 
value of kp is given a fixed value before the integral equation is to be solved. For values of 
kp up to some critical value of (A:i;') critical i the naive method successfully solves the integral 
equation within a finite number of iterations. 

In order to investigate the convergence we use a cobweb diagram, as shown in Fig. 4.4 in 
which we 'join-the-dots' between successive guesses and function evaluations. For clarity, the 
function f{M*,kp) is shown, as is the value of M*; the point at which these lines intersect 
represents the solution to the integral equation. In this way we can trace the convergence of 
the guesses from an initial guess of M* = to the value that satisfies the integral equation. 

If however, we wish to solve the integral equation for a value of kp > (A;_f) critical) the 
naive method fails to converge, and rather 'fiip-flops' between a pair of points. The cobweb 
diagram for this case is shown in Fig. 4.5, and we see that the solution is never reached. 



4.2. Steffensen's Method 



59 



> 

■ I— I 

1^ 



3 



10 

-1 



10 
10" 
10" 
10" 
10" 
10" 
10" 
10" 
10" 
10 



-10 



1 1 1 1 

5 iterations 
10 iterations 
' 15 iterations 
30 iterations 



0.0 



0.2 



0.4 



0.6 



0.8 



1.0 



m 
PI 

0) 

«4-H 

0) 







10 
-1 



10 
10" 
10" 
10" 
10" 
10" 
10" 
10" 
10" 
10 



-10 



5 iterations _ 
11 iterations 



0.0 



0.2 



0.4 



0.6 



0.8 



1.0 



Fig. 4.3: (Color Online) Errors for convergence of a naive method (top) and Steffensen's 
method (bottom) relative to the exact solution u*{x). The naive method converges in 20 
iterations, while Steffensen's method converges in 11 iterations. 



The significance of the value of (A:ir) critical is that the slope of the function at the fixed- 
point becomes too negative; in particular, /'(M*, (A;i7')criticai) < —1- Proving that this will 
cause this method to fail to converge is as easy as creating a cobweb diagram for any function 
with /'(fixed-point) < — 1. Even if we select a starting point particularly close to the actual 
fixed-point solution and perform this iterative procedure, if kp > (A:ir) critical the iterations 
diverge to a pair of points, as shown in Fig. 4.6. Although we have only shown the case of 
the naive substitution iterative procedure, we note that Newton's method suffers from the 
same limitation. 

If however, we use Aitken-improved Steffensen's method for solving the same integral 
equation with the same initial value of M* = 0, the procedure converges to the correct result, 
irrespective of the choice of kp. Fig. 4.7 shows an example of this convergence for a value 
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1 2 3 4 5 

M* [fm-^ 

Fig. 4.4: (Color Online) Cobweb diagram for a naive method of solving Eq. (4.10) for a single 
baryon, in which kp < (^F)criticai- The short-dashed line traces the guesses and function 
evaluations from the initial guess of M* = to the value that satisfies the integral equation. 

of kp > (kp) critical) the solution is not only found, but is done so in a small number 
of iterations. We shall therefore rely on Aitken-improved Steffensen's method to solve the 
integral equations in the work that follows. 

While we have shown an example of a single integral equation, the equations we will 
need to solve for our calculations will be a series of coupled integral equations, = 
/(M*, M*, kpp, kpj) and this complicates matters. To be precise, we would require an equiv- 
alent of the Jacobian, Eq. (4.9) but derived in the manner of Steffensen's method. For our 
purposes, we are able to use the linear version to solve two coupled integral equations; one 
for protons and one for neutrons. This may be due to the similarities between the values of 
M* for each of these. Further research will be undertaken to investigate the roles that the 
octet baryons take in these calculations. 
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Fig. 4.5: (Color Online) Cobweb diagram for a naive method of solving 
Eq. (4.10) for a single baryon, in which kp > (^F)criticai- In this case, the 
procedure does not converge to the a single solution, but alternates between a 
pair of points. 
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Fig. 4.6: (Color Online) Cobweb diagram for a nai've method of solving 
Eq. (4.10) for a single baryon, in which kp > (A;F)critical and the starting point 
is deliberately chosen to be near the exact solution. In this case, the procedure 
does not converge to the a single solution, but in fact diverges, then alternate 
between a pair of points. 
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Fig. 4.7: (Color Online) Cobweb diagram for a Aitken-improved Steffensen's method of solv- 
ing Eq. (4.10) for a single baryon, in which kp ^ (A;_p) critical) with initial guess M* = 0. This 
shows that even for cases in which the naive method (and Newton's method for that matter) 
would fail, this method is sufficient to solve the integral equation. 
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4.3 Infinite Matter 

Utilising the above techniques we are able to compute simulations of infinite matter. The 
entire procedure is coded in Fortran 90. The method by which we do this is as follows; 

1. Define all vacuum masses, coupling constants, and internal parameters. 

2. Provide suitable initial guesses for variables kp^^ and kp^^, and if applicable, {a). 

3. Select a value for the control parameter ptotai- We select a small non-zero value to begin 
with. When we are computing models involving the leptons, we also select the value of 
the total charge, which is always Qtotai = 0. 

4. Call a subroutine SolveConstraints which solves the constraint equations Eqs. (4.6)- 
(4.7) using Newton's Method. Eq. (4.8) can be excluded here even if not calculating 
Fock terms, but this self-consistency must then be solved independently. This subrou- 
tine takes initial guesses for kp^ and kp^ as well as as required tolerances as inputs, 
and produces the Fermi momentum values that satisfy the constraint equations above 
as outputs. 

5. As part of the evaluation of Eqs. (4.6)-(4.7) within SolveConstraints, the self-consistent 
effective masses must be calculated by calling another subroutine SolveEf f ectiveMasses 
This subroutine takes initial guesses for the effective masses and self-energies, as well 
as required tolerances as inputs, and produces the self-consistent effective masses that 
satisfy Eq. (4.8) as outputs. 

6. With the Fermi momenta and effective masses determined, we finally have all the in- 
formation required to evaluate quantities of interest, such as the energy density (e.g. 
Eq. (3.15)) and pressure (e.g. Eq. (3.16)) of the system. These values are evaluated 
and saved to a data-file. 

7. Finally, the control parameter ptotai is increased by some small step size 6p and the 
process repeated, where the initial guesses for the subroutines are now the final values 
from the previous iteration. 

This of course assumes that the couplings constants are known quantities. In each of the 
models we use we require that the properties of saturated nuclear matter are reproduced, 
as per Eq. (2.73). In order to find the appropriate values of g^a and gNui that reproduce 
the correct results we perform the above method with initial guesses for these values, and 
calculate the energy per baryon at saturation (as defined by the lower extrema of the data) 
and the Fermi momentum at which this occurs. By performing this several times for several 
choices of couplings, we can refine our choices for initial guesses in order to better reproduce 
the saturation properties. Of course, in practice we perform these steps using a Fortran 90 
subroutine which minimises a function of the saturation property variables for values of the 
couplings, such that the couplings found using this method produce saturation properties for 
any EOS that match the constraints with great precision. 

This process is suitable for solving all the models presented in this thesis; the model- 
specific details do not require changes to the method of solving the equations. This provides 
us with enough data to investigate fully the properties of infinite matter. To investigate 
compact stellar objects however, we require a further calculation. 



64 



4.4. Runge-Kutta Integration 



4.4 Runge-Kutta Integration 



In order to calculate the properties of compact stellar objects, we need to solve the Tolman- 
Oppenheimer-Volkoff Equations, as described in Section 2.11 and derived in Appendix A. 4. 
These equations provide a connection between the infinite matter equation of state and 
the properties of macroscopic stellar objects which can support the enormous gravitational 
pressures involved without collapsing. 

The relation between the state variables of the EOS; the energy density £ and the pressure 
P; and the the mass-radius relationship of a compact stellar object is determined by the 
pressure gradient 

dP _ G {P/c^ + £) {M{r) + 47rr3p/c2) 
Ik ~ r(r - 2GM(r)/c2) ' 

In order to make use of this, we will need to integrate this equation to find the pressure 
at some radius R' . We use fourth-order Runge-Kutta integration, which provides a useful 
scheme for integrating since it only requires a small number of function evaluations. 

Fourth-order Runge-Kutta integration is defined by the following iterative procedure; If 
the differential equation y' = G{x, y) is to be integrated to give 

y{h) = Vo+ f G{x,y)dx (4.19) 
Jo 

with initial values and intervals of 

yo = y(0), xo = 0, Xi = i5x, (4.20) 
where b = nJx, then we can define the Runge-Kutta steps as 

RKli = G{xi,yi), 

(OX ox \ 

RK^i = G (^Xi + y , + y RK2^ , 

RKAi = G{xi + 5x, yi + 6xRK3i), (4.21) 

where 6x is a small shift in x that will determine the resolution to which we wish to know 
y{x). The values RKli and RKAi are the slopes of tangents to the function y at the points Xi 
and Xi + 5x respectively, while RK2i and RK2>i are the slopes of y at the midpoint Xi + 5x/2, 
where the y value of the former is defined by the Euler method'^ using RKli, and the latter 
defined by the Euler method using RK2i. 

The next value of the function, located at Xj+i = (i + 1) 5x, is then defined using a 
modification of Simpson's Rule, as 

5x 

Vi+i = yi + -^ (RKli + 2 RK2i + 2 RK^i + RKA,) . (4.22) 



''Euler's method provides a very simple approximation, utilising the first two terms of a Taylor expansion 
of y to solve y' {t) = f{t,y{t)), thus y„+i = yn + hf{tn,yn). 
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By performing this procedure iteratively we can integrate a differential equation until after 
n steps we reach the upper integration limit of x„ = b. 

For our purposes, we select a central density of a stellar object /Occntrai as the control 
parameter. The EOS data calculated for a particular model then provides the energy density 
and pressure at the centre of the star. 

Assuming spherical symmetry, we require that the pressure gradient dP/dr at the centre 
of the star (r = 0) be zero, thus y'(0) = 0. By using Runge-Kutta integration, we can 
evaluate the pressure at the boundary of a sphere of matter with a radius r = R' . By relating 
this pressure to the EOS we can also find the energy density at this boundary, £{R'). 

By performing this integration for larger and larger values oi r = R' (since Eq. (4.18) 
describes a negative pressure gradient, this results in smaller and smaller values of pressure 
P) we eventually find some value r = R for which the integral of Eq. (4.18) becomes negative, 
i.e. Vi+i < 0; thus at a distance of approximately r = R (to within the resolution of 5x) the 
pressure is zero, and we define this point as the edge of the star, and hence the radius. 

Since we have also collected data for £{r), we can calculate the mass residing within some 
radius R' as 



and thus the total mass of the star is defined by M(R), where R is the radius of the star. With 
this information we can investigate the mass-radius relation for a star with a given central 
density. Repeating this calculation for several values of Pcentrai provides a locus of values 
as shall be shown, for example, in Section 5.1.2. We can also investigate other properties 
of a star as a function of internal radius by calculating them at various values of Xj = r in 



In order to combine the EOS for hadronic and quark matter (since we will be using a Gibbs 
transition, refer to Section 2.9) we require a method of calculating the phase transitions from 
hadronic matter to a mixed phase, and from a mixed phase to quark matter. Following the 
requirements described in Section 2.10 we calculate the EOS of hadronic matter with control 
parameter ptotah while at each density we use fin and fie as inputs to the quark matter EOS 
(for a given model of quark matter). At each density we can calculate the pressure for the 
hadronic phase Ph and the quark phase Pq and compare these values. At the point (if it 
exists) where Ph = Pq we consider the models to be in equilibrium. 

We then change the control parameter in our code from ptotal to x (parameterizing the 
quark fraction in the mixed phase, and acting as the order parameter between the phases) at 
which point we have x = 0, still purely hadronic matter. We can then increase x by Sx and 
calculate properties of the mixed phase such as the total baryon density and energy density 
with weightings of x and (1 — x) for the quark and hadronic phases respectively, with the 
condition that the chemical potentials are still equal for each phase, and the pressures are 
equal for each phase. 

We continue to increase x until we reach x = 1 at which point the EOS is described 
entirely by quark matter, and thus the mixed phase is in equilibrium with a pure quark 
phase. From this point on, we return to using ptotal as the control parameter and increasing 
this value by 6p up to some arbitrary density. 




(4.23) 



Eq. (4.20). 
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Results 



Within this chapter we will once again use the term 'configuration of a model' to indicate 
differences (such as types of particles included or neglected) within a particular model. Some 
of the most interesting and important numerical results are shown, though results for all 
calculations for all configurations will not be shown, due to the overwhelmingly large number 
of possibilities that exist. A summary table is provided in Section 5.5 which contains many 
numerical results of interest. 

5.1 QHD Equation of State 

Though QHD has been studied extensively, and many excellent summaries exist (e.g. Ref. [18]), 
we will present the results of our QHD calculations for the purpose of comparison in later 
sections, and to verify that our results do indeed reproduce the established results. 

5.1.1 QHD Infinite Matter 

To obtain numerical results, we solve the meson field equations, Eqs. (3.6)-(3.8), with the con- 
ditions of charge neutrality and fixed baryon density for various configurations of models. For 
configurations involving leptons, we include the condition of equivalence of chemical poten- 
tials given by Eq. (3.18). The energy per baryon given by Eq. (2.73) for various configurations 
of nucleonic QHD (in which QHD-I neglects contributions from p mesons, QHD-II includes 
these contributions, nuclear QHD contains equal proportions of protons and neutrons, and 
in each case we do not model leptons) are shown in Fig. 5.1 and we see that the saturation of 
nuclear matter occurs at the correct value kp (here we use the neutron Fermi momentum, 
as it is common to all three configurations; for the relation between Fermi momentum and 
density, refer to Section 2.4) corresponding to po as per Eq. (2.38) by construction via the 
use of appropriate couplings g^a and qnlu- 

The compression modulus for symmetric nuclear matter, as defined by the curvature at 
saturation, refer to Eq. (2.75) is found to be X = 525 MeV, which is in agreement with 
Ref. [18], but as stated in that reference, not with experiment. This point will be discussed 
further in Sec. 5.2.1. We further note that the pressure of nuclear matter calculated via 
Eq. (3.16) is zero at the saturation point, as predicted by Eq. (2.71). 

The case of 'symmetric' nuclear matter — in which we include protons, but not /3-equilibrium 
with leptons — is of course purely academic, since this is infinite matter, and we are therefore 
considering an infinite charge. The saturation that occurs in this case and the binding (a 
negative energy per baryon indicates a binding energy) is considered to resemble heavy finite 
nuclei. 

The low density EOS for nucleonic QHD is shown in Fig. 5.2 in which we see that the low 
density EOS is rather soft (does not approach the limit of P = £; the stiffest possible EOS in 
which the speed of sound equals that of light), and nearly configuration independent in that 
QHD-I produces similar results to that of nuclear QHD. At higher densities, the EOS for all 
configurations approach the limit. 
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Fig. 5.1: (Color Online) Energy per baryon for nucleonic QHD-I, QHD-II, and nuclear QHD 
which demonstrates saturation at the correct value of E as defined in Eq. (2.73), and value of 
kp (here we use the neutron Fermi momentum, as it is common to all three configurations) 
corresponding to pQ as per Eq. (2.38). The reproduction of this value occurs by construction 
via the use of appropriate couplings g^^ and qnu) as described in Section 4.3. 



The effective masses for the various configurations of nucleonic QHD are shown in Fig. 5.3, 
for which we note that the neutron matter curves are identical, as the p meson (a vector me- 
son) provides no contribution to the effective mass, which is a purely scalar effect. Although 
the form of the effective mass is linear in Fig. 3.1, in that case the effective mass is plotted 
against the scalar self-energy, which is shown in Fig 5.4, and we note that this quantity is 
non-linear, which determines the shape of the curve in Fig. 5.3. For the case of nuclear QHD 
we find values of the effective mass at saturation of the neutrons to be (M*/M„)sat = 0.56, 
corresponding to an effective mass of M* = 526 MeV. 

We can investigate the balance between the various meson fields by examining the self- 
energy contributions, which are defined explicitly for QHD in Eq. (2.42). These are shown in 
Fig. 5.4 for QHD-I and QHD-II both in /3-equilibrium (the self-energies for nuclear QHD are 
the same as those for QHD-I, since the p meson does not contribute to nuclear QHD as its 
contribution is proportional to the asymmetry between proton and neutron densities). From 
this figure it is clear that there exists an important balance between the scalar and vector 
interactions. Note the significance of the shape of the curve to that of the effective mass 
as shown in Fig. 5.3. Also note that in Hartree-level QHD, S^, = for all baryons. 

A log-log graph of the EOS for nucleonic QHD in which we show a Maxwell transition 
for the van der Waals style liquid-gas phase transition in QHD-I is shown in Fig. 5.5. The 
transition pressure is taken from Ref. [IS]. In this form of a transition (in contrast to that 
of Section 2.9) the pressure is constant between the phases (isobaric transition) and the 
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Fig. 5.2: (Color Online) Low density EOS for nucleonic QHD, as well as the causal limit of 
P = f in which the speed of sound in-medium equals the speed of light. Here we see that the 
low density EOS are rather similar, independent of which configuration of QHD is considered. 
The common feature however is that the low density EOS are rather soft, in that it does not 
yet approach the limit oiP = E. The EOS for all densities for these configurations are shown 
in Fig. 5.5. 



density of the combined EOS becomes disjoint. This example refers to a liquid-gas style 
phase transition, where later we will discuss phase transitions between phases with different 
degrees of freedom. The inclusion of the p meson removes this transition. 

For nuclear QHD, the species fractions (refer to Eq. (2.40)) are by definition equal for the 
protons and neutrons, i.e. Yp = Yn = 0.5. For QHD-I and QHD-H the neutrons provide the 
only baryonic contribution. If we include protons and leptons {i = e~, /x~) we can investigate 
the effects of /3-equilibrium matter with global charge neutrality. The species fractions Yi are 
shown in Fig. 5.6 for a configuration with (and without) contributions from p mesons. In 
each case, at low densities the system is composed of nearly entirely neutron matter, but we 
can see that at higher densities the contributions from protons become non-negligible. 

We do not show any results for the addition of hyperons to the QHD EOS. The QHD 
model as described in Section 3.1 is derived assuming a weak-field limit (in order to define the 
MFA; refer to Section 2.2) in order to define the Dirac equation for the baryons. In this limit, 
where the meson fields are assumed to be small in magnitude (weak) , we can ignore the effects 
of antibaryons, since it is not possible in this limit for the creation of particle-antiparticle 
pairs via vacuum fiuctuations^. The inclusion of hyperons causes the magnitudes of the fields 
— particular the scalar field, as can be seen by examining the form of this in Eq. (3.17) 



^This is a well known limitation of the Dirac equation, and has been discussed in detail, for example in the 
case of electron scattering from a potential barrier, resulting in the 'Klein Paradox' [lis]. 
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Fig. 5.3: (Color Online) Neutron effective masses for the various configurations of nucleonic 
QHD. Note that the curves for neutron matter (QHD-I and QHD-II) are identical; the addi- 
tion of the p meson contribution (a vector meson) does not contribute to the effective mass, 
as this is a purely scalar effect. Although the form of the effective mass is linear in Fig. 3.1, 
in that case the effective mass is plotted against the scalar self-energy, which is shown in 
Fig 5.4, and we note that this quantity is non-linear, which determines the shape of the curve 
here. 

in which we sum over all baryons — to become large, and thus the weak-field assumption 
becomes violated. We can observe that this is indeed the case by examining the effective 
masses. We observe that when the hyperons are included, the baryon effective masses for 
several of the baryons become negative, indicating that the scalar self-energy (defined by the 
mean-scalar field as per Eq. (3.10)) has become larger than the vacuum baryon mass. We 
consider this to be a breakdown of the model, and thus we do not perform calculations where 
this occurs. Calculations for QHD can be performed prior to this breakdown, but we elect to 
not present any of these results, as they represent a model which we consider to be inaccurate. 
We shall therefore consider it not possible to include hyperons into QHD in an interesting 
manner for our purposes. The issue of negative effective masses is however remedied in QMC 
as will be shown later. 
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Fig. 5.4: (Color Online) Self-energies for /3-equilibrium QHD-I and QHD-II. Nuclear QHD 
has the same self-energies as QHD-I since the p meson does not contribute; it is proportional 
to the asymmetry between proton and neutron densities. Here the density is shown up to 
1.0 fm~^ after which the curves become linear. Note the significance of the shape of the 
curve to that of the effective mass as shown in Fig. 5.3. Also note that in Hartree-level QHD, 
= for all baryons. 
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Fig. 5.5: (Color Online) Log-log EOS for nucleonic QHD-I and QHD-II. QHD-I produces 
features similar to that of a van der Waals EOS, which can be interpreted as a liquid-gas 
phase transition. The addition of the p meson contribution removes this transition. A 
Maxwell construction has been performed (using the transition pressure from Ref. [18]) to 
produce a constant pressure transition between the phases. Note that at high densities (large 
values of £) all the configurations approach the limit oi P = £. 
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Fig. 5.6: (Color Online) Species fractions Yi for nucleonic QHD in /3-equilibrium with leptons 
i = e~,/i~, including (left) and lacking (right) a contribution from the p meson. In both 
cases, at low densities the system is composed of nearly entirely neutron matter, but at higher 
densities the contribution from protons become non-negligible. 
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5.1.2 QHD Stars 

Having found the QHD EOS by evaluating the energy density, Eq. (3.15), and pressure, 
Eq. (3.16), we can calculate properties of stellar objects (which we shall refer to as 'stellar 
solutions') based on this EOS, using the TOV equation (refer to Section 2.11). The radius 
of a star is defined as the radius at which the pressure is zero and is calculated using a 
fourth-order Runge-Kutta integration method (refer to Section 4.4). 

The mass-radius relations for various configurations of nucleonic QHD are shown in 
Fig. 5.7 along with the 2a results^ of various experiments [69-71] — which detail the only 
measurements to date of the radii of neutron stars along with their masses, namely for the 
stars denoted by EXO 0748-676, 4U 1608-52, and EXO 1745-248~for comparison. The au- 
thors of Ref. [69] claim that their findings (corresponding to EXO 0748-676) rule out soft 
equations of state, though as noted in a response to that paper [72] and in this work, the 2a 
data admits a wide variety of EOS. The results of Ref. [69] should likely be discounted, since 
developments in the field have shown possible large errors with that experiment [73]. 




5 10 15 20 25 

R [km] 

Fig. 5.7: (Color Online) Mass-radius relation for various configurations of nucleonic QHD. 
Also shown are the theoretical limits due to General Relativity and causality (where G is the 
gravitational constant, c is the speed of light), as well as the 2a error bounds from various 
experiments [69-71] which illustrate how loosely the requirements of an EOS are currently 
defined. Note that the low-density shape of the curve changes with the addition of the p 
meson in QHD-II. 

^The shaded areas for EXO 1745-248 and 4U 1608-52 represent a conservative reproduction of the 2a data 
from the relevant reference. The shaded area for EXO 0748-676 represents the central data point plus error 
bars expanded to a rectangular area. 
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We note that the shape of the low-central-density relation changes with the addition of 
the p meson producing QHD-II, which should be expected, given the differences between the 
QHD-I and QHD-II EOS as shown in Fig. 5.5. 

We show (as an example, for QHD-I) the relation between pressure and internal radius 
for a variety of stars with different central densities in Fig. 5.8. We note the trend that a star 
with a larger central density has a greater central pressure, and a smaller total radius. Each 
of the curves in that figure correspond to the high-central-density positive-gradient section 
of the relevant curve in Fig. 5.7. 

If we now include the leptons to our calculations, we can observe the effect that these have 
on the mass-radius relations, as shown in Fig. 5.9 for /3-equilibrium nucleonic matter. While 
the QHD-I and /? QHD-I configurations may not be dramatically different, the similarities 
of the stellar solutions for these configurations demonstrates the difficulties that exist in 
determining the content of stars from these two parameters alone. 




R (internal) [km] 

Fig. 5.8: (Log) pressure vs. internal radius for stars with various central densities for QHD-I. 
Note that the stars with a larger central density have a larger central pressure, and a smaller 
total radius. The data for Pccntrai = 0.5 fm~^ corresponds to a star with M = 2.84 Mq 
and R = 13.96 km, while the data for Pcentrai = 1-75 fm~^ corresponds to a star with 
M = 2.39 Mq and R = 10.9 km, both of which correspond to the high-central-density 
positive-gradient section of the relevant curve in Fig. 5.7. 
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Fig. 5.9: (Color Online) Mass-radius relations for QHD-I and QHD-II in /3-equilibrium, with 
the neutrons-only relations (as shown in Fig. 5.7) for comparison. Note that the effect of 
including protons and leptons does not have such a pronounced effect on QHD-I, since in 
that case the protons do not provide such a large contribution, as can be seen in Fig. 5.6. 
The similarities between the stellar solutions for each configuration emphasises the difficulty 
of using the mass and radius alone to constrain a model. 
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5.2 QMC Equation of State 

By including the effect of the quark content of baryons we can investigate the impact on 
various quantities, using the QMC model as described in Sec. 3.2. We will see that this 
has profound consequences for the EOS. Many of the results presented in this section are 
published by the author [42] as a unique investigation of the octet QMC model. 



5.2.1 QMC Infinite Matter 

The energy per baryon curves for QMC are shown in Fig. 5.10 where, for the sake of com- 
parison to Fig. 5.1 we have used the same notation for distinguishing the configurations, thus 
QMC-I neglects the contribution of the p meson, QMC-II includes it, and QMC nuclear re- 
stricts the species fractions of protons and neutrons to be equal. We neglect the contributions 
of leptons for these simple configurations. The couplings of the baryons to the mesons are 
found such that once again the appropriate saturation properties are reproduced. 

An interesting feature to note for Fig. 5.10 is that the curvature for nuclear QMC at 
saturation is less than that of the same curve in Fig. 5.1 indicating that the compression 
modulus for QMC takes a smaller value. This is of particular interest, since the value of 
the compression modulus for QHD is known to be too large [18], which is typical of models 
that neglect quark level interactions. For QMC we find a significant improvement in the 
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Fig. 5.10: (Color Online) Energy per baryon for nucleonic QMC-I, QMC-II, and nuclear QMC 
which demonstrates saturation at the correct value of E as defined in Eq. (2.73), and correct 
value ol kp. For comparison to Fig. 5.1 we use the same notation to differentiate between 
the configurations. The reproduction of the saturation value occurs by construction via the 
use of appropriate couplings and g^^j which have been re-calculated for QMC. 
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compression modulus; K = 280 MeV which hes at the upper end of the experimental range. 
The nucleon effective mass at saturation for QMC is found to be (M*)sat = 735 MeV, 
producing (M*/M)sat = 0.78. The effective masses for the various configurations of nucleonic 
QMC are shown in Fig. 5.11 and we see that these effective masses do indeed remain positive, 
as opposed to the negative effective masses we encountered in QHD. 

The most interesting aspects of QMC are found when we include the remainder of the 
octet of baryons — the hyperons — which was not possible in QHD due to the Klein Paradox. 
As discussed in Sec. 3.2, the baryon effective masses in QMC do not become negative at any 
density, and thus we do not encounter the same issues as we do when trying to calculate octet 
QHD. The species fractions for octet QMC are shown in Fig. 5.12, where we note that the 
A species fraction is significantly larger than that of the S. The investigations by Rikovska- 
Stone et al. [41] (in which alternative Fock terms were introduced to the QMC model) lead 
us to expect that the S would disappear entirely from the system if we were to include Fock 
terms for the hyperons. 

The species fractions in Fig. 5.12 are particularly interesting when compared to the (3- 
equilibrium nucleonic configurations (Fig. 5.13) as the leptonic contributions decrease at 
higher densities, in contrast to the plateau observed in Fig. 5.13. This is a result of the 
chemical potential equilibria which provide preference to hadronic charged states with large 
masses, rather than those with smaller masses. In this case, the system is more stable with 
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Fig. 5.11: (Color Online) Effective neutron masses for the various configurations of nucleonic 
QMC in which the effective masses have a quadratic form. Note the subtle differences in 
shape between these curves and those of Fig. 5.3, particularly the value of the effective mass 
at saturation. 
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Fig. 5.12: Species fractions Yi for octet QMC, in which the p meson contribution is included. 
The density here is the total baryon density, as per Eq. (2.38). Note that in this case, all of 
the octet baryons contribute at some density, and that the species fractions of S hyperons are 
suppressed — in particular that of the S~ — compared to the other baryons. When compared 
to the nucleon-only data (as in Fig. 5.13) we observe that in this case the lepton densities do 
not plateau, but rather decrease once the hyperons appear, as explained in the text. Note that 
the A hyperon is the first hyperon to appear, and that at high-densities the relative species 
fractions plateau for the nucleons. A, and H hyperons. The relative baryon proportions at 
high density are ordered by isospin, in that pA > Pa^ ~ Ph > Pe- The parameters used here 
are shown in Table 5.2. 
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Fig. 5.13: (Color Online) Species fractions Yi for nucleonic QMC in /3-equilibrium with leptons 
i = e~ , n~ , including (left) and lacking (right) a contribution from the p meson. These results 
are very similar to those of Fig. 5.6, though the differences are most pronounced when lacking 
the p meson (right). 
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the negative charge (to counter that of the protons) being provided by the H~ hyperon rather 
than the leptons. 

One may reasonably ask then why the T,~ meson does not play the role of balancing pos- 
itive charge, as it is the lightest baryon with a negative electric charge. The relation between 
the chemical potentials in Eq. (2.46); /i^- = f^E- involves the full interacting chemical poten- 
tials of Eq. (2.44) which includes scalar meson terms that affect the mass, and vector meson 
terms that effect the energy. It is a balance of these effects which dictates which particle 
will have the greater Fermi momentum (hence, density) when considering Eq. (2.46) above. 
For the case shown in Fig. 5.12 the balance dictates that the H hyperons will have a greater 
species fraction at the densities shown. 

If however we exclude the contribution of the p meson, we can observe the effects that 
are attributed to this meson. The species fractions for octet QMC in which we neglect the p 
meson are shown in Fig. 5.14, and we note that the distribution is remarkably different. First, 
we note that the S~ hyperon is now the first to appear, in stark contrast to Fig 5.12 in which 
the S hyperons were largely suppressed. The effect on the chemical potentials (in particular 
the vector potential terms) due to neglecting the p meson is such that this is now possible. 
We also note that the high-density plateau of Fig. 5.12 is less evident, and the baryons are 
not sorted by isospin as distinctly. The lack of suppression for the S~ is made clear in this 
case, and indicates the strong link between this meson and baryon, which is reasonable given 
that the p meson couples to isospin, and the S baryons have the largest magnitude isospin 
(/e = 1). 




Density p [fm ^] 



Fig. 5.14: Species fractions Yi for octet QMC, in which the p meson contribution is neglected, 
but otherwise using the same parameters as used in Fig. 5.12. We note that at high-densities, 
the baryons are no longer sorted by isospin as distinctly, and that the S~ hyperon is now 
the first to appear. Furthermore, the hyperon threshold now occurs at a higher density than 
that of Fig. 5.12 where the A was the first to appear. In this case, the S~ species fraction is 
not suppressed at high density, indicating that the p meson plays a vital role in this. 
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5.2.2 QMC Stars 

By solving the TOV equation, we can once again investigate 'stellar solutions'. The mass- 
radius relations for nucleonic QMC EOS are shown in Fig. 5.15. The maximum masses 
for the QMC-I and nuclear QMC EOS are much smaller than their counterparts in QHD, 
since the QMC EOS is much softer than that of QHD as evidenced by the smaller value of 
the compression modulus K in QMC. The softening is due to the effective inclusion of the 
response to the scalar field via the scalar polarizability, as discussed in Section 3.2. Inclusion 
of the p meson appears to undo much of this softening. We do note however that the QMC-I 
curves still lie within the 2a bounds. 

If we now introduce hyperons to the QMC EOS, as shown in Fig. 5.16 we see a further 
softening of the EOS. In this case, the addition of the p meson does not stiffen the EOS back 
to the level of the nucleonic EOS, since the p contribution is much smaller in the hyperonic 
case. As shown in Fig. 5.12, the asymmetry between the components of the isodoublets and 
those of the isotriplet is small, and the p contribution is proportional to this asymmetry (refer 
to Eq. (3.8)). This is an important distinction between the nucleonic and hyperonic EOS. 

While the mass-radius relations for these configurations may be incapable of explaining 
the given observed data of [69], it is clearly a great step forward in that we may now model 
the effects of hyperons, which was not possible in QHD. 




R [km] 

Fig. 5.15: (Color Online) Mass-radius relations for various configurations of nucleonic QMC, 
with additional features the same as Fig. 5.7. Note that the maximum masses of QMC-I 
and nuclear QMC stars are much lower than for QHD-I and nuclear QHD, indicating that 
the QMC EOS is much softer that that of QHD as evidenced by the smaller value of the 
compression modulus K in QMC. 
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Fig. 5.16: (Color Online) Mass-radius relations for octet QMC-I (as per Fig. 5.14) and QMC- 
II (as per Fig. 5.12), with additional features the same as Fig. 5.7. Note that the maximum 
masses of octet QMC-I and octet QMC-II stars are lower than nucleonic QHD-I and nucleonic 
QHD-II stars, indicating that the octet QMC EOS is softer that that of the nucleonic QHD. 



If we consider the case of a star with a central density of Pccntrai = 1-2 fm~ (an arbitrary 
choice) as shown in Fig. 5.17 for the case of hyperonic QMC, we observe that the outermost 
3 km of this 11 km star contains only nucleonic matter (in /3-equilibrium) . The core of this 
star contains roughly equal proportions of nucleons, A and H hyperons, with a notable lack 
of S hyperons. This is the simplest case for which we produce a stellar object with hyperons, 
and comprises our simplest method for including a strangeness degree of freedom. 
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Fig. 5.17: Species fractions for octet QMC in /3-equilibrium as a function of stellar radius for 
a stellar solution with a central density of /Ccentrai = 1-2 fm~^. The parameters used here are 
the same as those used to produce Fig. 5.12. We note that the outermost 3 km of this star 
contains only nucleons and leptons, with hyperon effects occurring much deeper. 
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5.3 Hybrid Equation of State 

A focus of this work has been to calculate the properties of hybrid stars; involving contri- 
butions from both baryons and deconfined quarks in a statistical mechanics method. The 
process by which we obtain these results has been described in Section 2.9. The results pre- 
sented in this section have been published by the author [42] as a novel extension of the octet 
QMC model. 

5.3.1 Hybrid Infinite Matter 

In Section 3.3 we discussed the MIT bag model for quark matter, in which we model three 
deconfined free quarks in a Fermi gas possessing constant masses consistent with current 
quark phenomenology. It is this model that we shall refer to (unless stated otherwise) for 
hybrid models involving quark matter. The hadronic model we will use in this section is the 
octet QMC model, unless otherwise specified. 

The conditions for a Glendenning-style mixed phase as discussed in Section 2.10 require 
that for a given pair of Hn and fie (common to the hadronic and quark phases) at any value 
of the mixing parameter x, the quark density is greater than the hadronic density. This 
condition ensures that the total baryon density increases monotonically within the range 
Pqp > Pmp > PHP, as can be seen in Eq. (2.92). An example of this is illustrated in Fig. 5.18 
for a mixed phase of octet QMC and three-flavor quark matter modelled with the MIT bag 
model. 

Recall that the method for evaluating the properties of a mixed phase involves calculating 
the neutron and electron chemical potentials of the hadronic phase and using these as inputs 
to the quark phase calculations, providing the quark Fermi momenta via Eq. (3.30). In this 
case, with the dynamic quark masses of NJL and no vector potentials (self-energies) due to 
non-interacting quarks, the calculated quark Fermi momentum corresponds a quark density 
that is lower than the hadronic density, and as a result there are no configurations for a 
mixed phase in which the proportion of quarks increases while at the same time the total 
baryon density increases. It may be possible that with smaller constituent quark masses at 
low density, the Fermi momenta would provide sufficiently high quark densities, but we feel 
that it would be unphysical to use smaller constituent quark masses. This result implies that 
— at least for the models we have investigated — dynamical chiral symmetry breaking (in the 
production of constituent quark masses at low density) prevents a phase transition from a 
hadronic phase to a mixed phase involving quarks. 

We do note, however, that if we restrict consideration to nucleons only within the QMC 
model (with the same parameters'^ as octet QMC), and represent quark matter with the NJL 
model, we do in fact find a possible mixed phase. More surprisingly, the density at which a 
phase transition to the mixed phase occurs for this combination is significantly larger than 
the case where hyperons are present (in which case we model a transition to MIT bag model 
quark matter). An example of this is shown in Fig. 5.19, the parameters for which can be 
found in Table 5.3. This produces a mixed phase at about 4pQ {p = 0.64 fm~^) and a pure 
quark matter phase above about lOpo {p = 1-67 fm~'^). In this case the u quark first appears 

■^The couplings of mesons to baryons are defined by tlie energy per baryon, wliicli saturates at a lower 
density than the hyperon threshold. We can therefore safely assume that the inclusion or removal of hyperons 
plays no part in defining the couplings of mesons to the nucleons. 
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Fig. 5.18: (Color online) Densities in the mixed phase for octet QMC mixed with three- 
flavor quark matter modelled with the MIT bag model. Note that at all values of x (the 
mixing parameter according to Eq. (2.92) which also defines the total density), the equivalent 
baryon density of quarks is greater than the hadronic baryon density, allowing the total baryon 
density to increase monotonically with increasing x- Clearly, at x = the density corresponds 
entirely to the hadronic phase, and at x = 1 the density corresponds entirely to the quark 
phase. Between these points, the density is that of the mixed phase. 



at a higher density than the d or s quarks due to its positive charge which can only be 
balanced by the other two quarks since the leptons provide a decreasing contribution as per 
Eq. (2.47), since the proton and neutron densities become more and more similar. Although 
this example does show a phase transition, the omission of hyperons is certainly unrealistic. 
This does however illustrate the importance and significance of including hyperons, in that 
their inclusion alters the chemical potentials which satisfy the equilibrium conditions in such 
a way that the mixed phase is no longer produced. 

Using the MIT bag model, we can calculate the hybrid EOS, which is shown in Fig. 5.20. 
We can see that at low densities (low energies) all of the configurations of EOS are fairly soft. 
At the lowest densities, each of the configurations of EOS are approximately equal; at this 
point they all represent nucleons in /3-equilibrium. 

For each of the configurations in which we find a phase transition from baryonic matter 
to quark matter, the EOS consists of negatively charged quark matter, positively charged 
hadronic matter, and a small proportion of leptons, to produce globally charge- neutral matter. 
The proportions of hadronic, leptonic and quark matter throughout the mixed phase (for 
example, during a transition from octet QMC matter to three-fiavor quark matter modelled 
with the MIT bag model) are displayed in Fig. 5.21 in which we note that the quarks are able 
to satisfy charge neutrality without lepton contributions, in contrast to the cases of nucleonic 
/3-equilibrium in which the lepton contributions remain stable at increasing densities. A 
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Fig. 5.19: Species fractions for a phase transition from nucleonic QMC matter to three-flavor 
quark matter modelled with NJL. Note that in this unphysical (we see no physical reason 
for which hyperons should be neglected) case, a phase transition is possible, and occurs at a 
value of p = 0.64 fm~^. In this case the u quark first appears at a higher density than the 
d or s quarks due to its positive charge which can only be balanced by the other two quarks 
since the leptons provide a decreasing contribution as per Eq. (2.47). 

summary of the results of interest is given in Table 5.3. 

Results of calculations for larger quark masses are not shown, as they require a much 
lower bag energy density to satisfy the pressure equilibrium conditions. For constituent quark 
masses, we find that no phase transition is possible for any value of the bag energy density, 
as the quark pressure does not rise sufficiently fast to overcome the hadronic pressure. This 
is merely because the masses of the quarks do not allow a sufficiently large Fermi momentum 
at a given chemical potential, according to Eq. (3.30). 

When we calculate the EOS including a mixed phase and subsequent pure quark phase, 
we find that small changes in the parameters can sometimes lead to very significant changes. 
In particular, the bag energy density B, and the quark masses in the MIT bag model have 
the ability to both move the phase transition points, and to vary the constituents of the 
mixed phase. We have investigated the range of parameters which yield a transition to a 
mixed phase and these are summarized in Table 5.3. For illustrative purposes we show an 
example of species fractions for a reasonable set of parameters (5^/^ = 180 MeV and mu,d,s = 
3, 7, 95 MeV) in Fig. 5.22. Note that in this case the A hyperon enters the mixed phase briefiy 
(and at a low species fraction). Note that the transition density of /?mp ~ 0.22 fm~^ produced 
by the combination of the octet QMC and MIT bag models in that case seems unlikely to be 
physical as it implies the presence of deconfined quarks at densities less than 2po , which would 
contradict the results of searches for such entities. A similar transition from nucleonic QMC 
matter to three-fiavor quark matter modelled with the MIT bag model (Fig. 5.19) produces 
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Fig. 5.20: (Color online) Equation of State for; nucleonic 'N' matter modelled with QMC; 
nucleonic matter where a phase transition to NJL modelled quark matter is permitted; bary- 
onic 'N+Y' matter modelled with octet QMC; and baryonic matter where a phase transition 
to MIT bag modelled quark matter is permitted. The line P = £ represents the causal limit, 
^^sound = c. The bends in these curves indicate a change in the composition of the EOS, such 
as the creation of hyperons or a transition to a mixed or quark phase. Note that at very low 
energies (densities) the curves are identical, where only nucleonic matter in /3-equilibrium is 
present. 



results almost identical to those of Fig. 5.22, except of course that in that case there is no 
contribution from the A hyperon. We note the significant differences in threshold densities 
between these two cases. 

With small changes to parameters, such as those used to produce Fig. 5.23 in which the 
bag energy density is given a slightly higher value from that used in Fig. 5.22 [B^^'^ increased 
from 180 MeV to 195 MeV, but the quark masses remain the same), it becomes possible for 
the H hyperons to also enter the mixed phase, albeit in that case with small species fractions, 
Yt.^Y'E. < 0.02. This highlights the need for strict tolerances for 'known' values of model 
parameters. 
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Fig. 5.21: (Color online) Charge-densities (in units of the proton charge per fm^) in the 
mixed phase for a transition from octet QMC to three-flavor quark matter modelled with the 
MIT bag model. The net charge of the mixed phase is exactly zero, as per the constraining 
equation of Eq. (2.95). Note that following the mixed phase, the quarks are able to satisfy 
charge neutrality with no leptons, in contrast to the nucleonic cases of /3-equilibrium in which 
the lepton contribution remains stable, x is the mixing parameter within the mixed phase 
according to Eq. (2.92). 
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5.3. Hybrid Equation of State 




Fig. 5.22: Species fractions Yi for octet QMC (the same as in Fig. 5.12) but where now 
we ahow the phase transition to a mixed phase involving quark matter modelled with 
the MIT bag model. Note that the A hyperon is the only hyperon to appear in the mixed 
phase, and does so at a much higher density than the configuration where the transition 
to a mixed phase is forbidden. A similar transition from nucleonic QMC matter to 
three-flavor quark matter modelled with the MIT bag model (as shown in Fig. 5.19) 
produces results almost identical to these, except of course that in that case there is 
no contribution from the A hyperon. We note the significant differences in threshold 
densities when comparing these results. 




Fig. 5.23: Species fractions, Yi, for octet QMC (the same as in Fig. 5.22 but now where 
the bag energy density has been increased to B^^^ = 195 MeV). Note that now the 
appearance of hyperons occurs at a smaller density than in the case of Fig. 5.22, the 
transition to a mixed phase occurs at a slightly larger density, and that now H hyperons 
are present in the mixed phase. 
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5.3.2 Hybrid Stars 

The stellar solutions for octet QMC hybrid stars (in which a phase transition to a mixed 
phase is allowed) are shown in Fig. 5.24 along with the baryon-only stellar solutions for 
comparison. We note that the solutions for configurations involving a phase transition are 
identical to those where a transition is neglected, up to some value of the central density at 
which point the solutions diverge. This indicates the lowest mass stars that contain quark 
matter. 

We note that the stellar masses for configurations involving a phase transition to quark 
matter are lower than those in which the transition is neglected, for the same value of central 
density. This is due to the softening of the EOS by the introduction of quarks (see Section 2.8). 
Overall, the stellar masses for these configurations are similar to observed neutron star masses, 
though notably lower than the masses of the most massive observed neutron stars. This could 
be attributed to an over-softening of the EOS. 

If we consider a star with a central density of /Occntrai — 

1.2 fm ^ for octet QMC with a 
phase transition to quark matter modelled with the MIT bag model (using the same param- 
eters as used in Fig. 5.22) we observe the prediction of 3.5 km of quark core for a star with a 
radius of 10 km, with only a very small contribution from the hyperons, via A, as shown in 
Fig. 5.25. 

If we increase the value of the bag energy density from i?^/^ = 180 MeV as used in 
Fig. 5.25 to B^/^ = 195 MeV we find that the radius of a star with a central density of 
Pcentrai = 1-2 fm~'^ increases by around 1 km, but the baryonic contribution (in the mixed 
phase) is more prolific, as we now find contributions from A and H hyperons in the core of the 
star, along with quark contributions. The effect of this increase in the bag energy can be seen 
by comparing Fig. 5.25 with Fig. 5.26, where the bag energy densities are S^/^ = 180 MeV 
for the former, and B^^^ = 195 MeV for the latter, and this is the only change that has been 
made. 
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5.3. Hybrid Equation of State 
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Fig. 5.24: (Color Online) Mass-radius relations for octet QMC in /3-equilibrium (neglecting 
the p meson, upper; and including it, lower) for baryonic and hybrid stars. Also shown are 
the additional features described in Fig. 5.7. Note that the stellar solutions for configurations 
involving a phase transition are identical to those where the phase transition is neglected, up 
to some value of the central density at which the solutions diverge. This point corresponds 
to the lowest density stars which contain quark matter. The stars containing quark matter 
generally have a lower mass than those which only contain baryonic matter due to a softening 
of the EOS. 
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Fig. 5.25: Species fractions for octet QMC with a phase transition to three-flavor quark 
matter modelled with the MIT bag model, as a function of stellar radius for a stellar 
solution with a central density of /Ocentrai = 1-2 fm~^. The parameters used here are 
the same as those used to produce Fig. 5.22. Note that in this case one finds pure 
deconfined three-flavor quark matter at the core (all of some 3.5 km) of this star, and a 
small proportion of A in the mixed phase. 
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Fig. 5.26: Species fractions for the interior of a star with central density Pccntrai = 
1.2 fm~'^ where the bag energy density is given a slightly higher value from that used in 
Fig. 5.25 (increased from B^^^ = 180 MeV to 195 MeV), but the quark masses remain 
the same. 
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5.4. Hartree-Fock QHD Equation of State 



5.4 Hartree-Fock QHD Equation of State 

As discussed in Section 3.5, we can further extend our models by including Fock terms. For 
comparison purposes, we show the effects of including only the Fock contribution to the 
self-energy, as well as the effects of the full Hartree-Fock calculations (including the medium 
polarization energy), as compared to the Hartree calculations presented earlier. 

The results presented in this section extend the calculations of Ref. [18] in which the 
properties of Hartree-Fock nuclear matter were calculated without /3-equilibrium. To the 
extent of the author's knowledge, the following calculations have not been performed or 
published elsewhere, and as such are further novel calculations. 

For simplicity, we present the extension from Hartree to Hartree-Fock for QHD, though 
similar calculations can be performed for extending QMC similarly. We leave this as work 
for the future. 

5.4.1 Hartree-Fock QHD Infinite Matter 

In Fig. 5.27 we show the saturation curves for nuclear QHD, in which we observe that the 
saturation properties for each configuration are reproduced accurately, in order to calculate 
the couplings to be used for each configuration. 

We also note from that figure that the curvature for each configuration is different at 
saturation, implying differences between the compression modulii for each configuration, with 
Hartree QHD having the largest value for K, and the full Hartree-Fock (Fock2) having the 
smallest. For comparison, the values are shown in Table 5.2. The softening of this EOS brings 
the compression modulus for Hartree-Fock QHD closer to the experimental range, though we 
note that this may further over-soften the EOS since the maximum masses in Hartree QMC 
appear too small already. 

We can examine the effective masses for Hartree-Fock QHD for the case of only including 
the self-energy Fock terms, and compare these with those for Hartree QHD as shown in 
Fig. 5.3. The comparison is shown in Fig. 5.28 and we note that the differences are subtle 
at this stage. The Hartree-Fock data corresponds to two configurations; Fockl denotes the 
inclusion of the self-energy Fock contribution, while Fock2 denotes the additional inclusion of 
the medium polarization energy Fock contribution. In each case, the effective mass is plotted 
for a limited range of Fermi momenta corresponding to densities up to ~ 0.5 fm~^. 
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Fig. 5.27: (Color Online) EOS for nuclear QHD matter calculated for three configurations; at 
Hartree level, including only Fock self-energy terms, and additionally including the medium 
polarization energy contribution. Couplings for each configuration are found such that sat- 
uration properties are reproduced. The values of these couplings can be found in Table 5.2. 
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5.4. Hartree-Fock QHD Equation of State 
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Fig. 5.28: (Color Online) Nucleon effective masses for nuclear QHD and /3-equilibrium nucleon 
QHD-I including Fock contributions to the self-energies (Fockl) and including this term as 
well as the medium polarization energy contribution (Fock2). The Hartree-Fock contributions 
are plotted for a limited range of Fermi momenta corresponding to densities up to ~ 0.5 fm^'^. 
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5.4.2 Hartree-Fock QHD Stars 

To investigate the effects that the addition of the various Fock terms have on stellar solutions 
for these EOS, we can solve the TOV equations in each case. Fig. 5.29 shows the mass-radius 
relations for nucleonic QHD-I with the additional self-energy Fock terms (Fockl) included, 
in which the couplings used provide a fit to the saturation properties of nuclear matter. For 
comparison purposes, the same stellar solutions with Fock terms neglected (as per Fig. 5.7) 
are shown, and we note that the maximum mass is lower when the Fock terms are included, 
indicating a softening of the EOS. 

Since the compression modulus of the full Hartree-Fock configuration (Fock2) for QHD 
is smaller than than of the Fockl configuration, we expect that maximum mass calculated 
with the Fock2 configuration would be lower than those in Fig. 5.29. 
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Fig. 5.29: (Color Online) Mass-radius relations for Hartree-Fock QHD-I, in which we have 
included the self-energy Fock terms. Also shown for comparison are the mass-radius relations 
of Hartree QHD-I as per Fig. 5.7. We observe that the inclusion of this Fock term softens 
the EOS, resulting in a reduced maximum stellar mass. 



96 
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5.5 Summary Tables 

The vacuum masses of the baryons, leptons, and mesons used in the calculations for this work 
are shown in Table 5.1, as per Ref. [17]. 

The couplings used in each configuration of each model are found such that the nuclear 
matter configuration reproduces saturation properties (as per Section 2.8). The couplings 
used for all of the individual configurations of QHD, QMC, and Hartree-Fock QHD are 
shown in Table 5.2. In the case of Hartree-Fock QHD we show the couplings that reproduce 
saturation properties for the case of only including the self-energy Fock contribution (denoted 
Hartree-Fockl), and for the case of including that contribution and additionally the medium- 
polarization Fock contribution (denoted Hartree-Fock2). The couplings are given to three 
decimal places to allow comparison between investigations, though more significant figures 
would be required to exactly reproduce these results. 

Table 5.3 contains comparison values for each of the configurations of each model; The 
input parameters used to calculate each configuration, such as species modelled, and — where 
applicable — the quark masses and bag energy density; as well as the calculated values for 
phase transition densities, compression modulii, and effective mass at saturation. 

Table 5.1: The vacuum (physical) baryon, lepton, and meson masses (in units of MeV) as 
used throughout this work [ ! 7] . 



Mp 
938.27 


Mn 

939.57 


Ma 
1115.68 


Ms- M^o 
1197.45 1192.64 


Ms+ 
1189.37 


M=- Mso 
1321.31 1314.83 




ma- 
550.0 




782.6 




nip 
775.8 






mg- 
0.51 




m^- 
105.66 





Table 5.2: Baryon-meson couplings for all configurations used in this work, found such that 
the saturation properties of nuclear matter (as described in Section 2.8) are reproduced. Also 
shown here are the resultant compression modulii for each configuration (we only need present 
the isospin symmetric value, as per the definition of K. 



Configuration 


gNa 


QNu) 


9p 


K (MeV) 


Hartree QHD 


10.644 


13.179 


6.976 


525 


Hartree QMC 


8.268 


8.417 


4.167 


281 


Hartree-Fockl QHD 


10.001 


11.819 




456 


Hartree-Fock2 QHD 


11.289 


14.170 




366 



Table 5.3: Table of species content (A^ = nucleons, Y = hyperons, i = leptons, q = quarks); inputs {B^^^, ^q)] and results for octet 
QMC and quark models presented in this thesis, py, pup and pqp represent the density at which hyperons first appear (A is the 
first hyperon to enter in these calculations in all but one configuration; octet QMC-I (Fig. 5.14 where the p meson is neglected, in 
which case the S~ is the first hyperon to enter); the density at which the mixed phase begins; and the density at which the quark 
phase begins, respectively. Figures for selected parameter sets are referenced in the final column. Dynamic NJL quark masses are 
determined by Eqs. (3.39)-(3.40). 



Particles: 








B^/^ (MeV) 


{mu,md,ms} (MeV) 


PY (fm-3) 


PMP (fm ^) 


PQP (fm ^) 


Fi^ 


fure: 


N, Y, e, 




UJ, 


p) 






0.33 






Fi^ 


5. 5.12 


N, Y, i, 




u 


) 






0.39 






Fij 


[. 5.14 


N, Y, i, q. 




UJ, 


p) 


180 


{3, 7, 95} 


0.54 


0.22 


0.95 


Fij 


[. 5.22 


N, Y, e, q. 




UJ 


) 


180 


{3, 7, 95} 


0.51 


0.36 


0.95 






N, Y, i, q. 




UJ, 


p) 


195 


{3, 7, 95} 


0.33 


0.35 


1.46 


Fi^ 


^ 5.23 


N, Y, i, q. 


{<r, 


UJ, 


p) 


170 


{30, 70, 150} 


0.55 


0.20 


0.87 






N, Y, i, q. 


{<r, 


UJ, 


p) 


175 


{100, 100, 150} 


0.41 


0.28 


1.41 






N, e, q. 




UJ, 


p) 


180 


Dynamic (NJL) 




0.64 


1.67 


Fi^ 


^ 5.19 
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Conclusions 



We have derived and performed calculations for various Quantum Hadrodynamics (QHD) 
equation of state (EOS) configurations (preliminarily at Hartree level) as a foundation for 
further work, both reproducing the results of a well-known reference, viz Ref. [18] and ex- 
tending the model to include /3-equilibrium between baryons and leptons. In this model the 
baryons are treated as effective degrees of freedom and interact with a mean- field of mesons. 
We have shown that this model is not capable of modelling hyperons due to a breakdown 
in the assumptions — in particular that the meson potentials are small, an assumption made 
in order to require the Dirac equation to be defined for this model — and thus we have not 
pursued this model any further, apart from calculating the stellar solutions; the mass-radius 
relations. 

We have investigated the latest (to date) manifestation of the quark-meson coupling 
(QMC) model to produce an EOS for nucleonic matter as well as investigations of the inclusion 
of hyperons to this model, in /3-equilibrium with leptons. In this model we include the 
self-consistent response of the internal quark degrees of freedom to the applied scalar field, 
allowing the hyperons to be modelled without violating the above assumptions of the Dirac 
equation that prevent us from calculating hyperon effects in QHD, as evidenced by the fact 
that the effective baryon masses in QMC — which now include a quadratic scalar-field term 
— remain positive at all densities. We do note however that other models for hadronic 
matter exist in which hyperon degerees of freedom become accessible, and in which case 
a transition from hadronic matter to quark matter becomes possible, e.g. Ref [74]. We 
have calculated the EOS for various configurations of the models described in this thesis with 
various values for variable quantities, and calculated constituents of infinite matter for several 
of these configurations, providing information on the density fractions of various particles in 
equilibrium. We have investigated the effects that changing various parameters of this model 
have on these density fractions of particles, as well as the effects of restricting the types of 
particles modelled with QMC. We note that qualitatively, the relative proportions of the 
hyperons are in good agreement with other treatments (for example, Ref. [75]) and that 
we observe expected phenomenology, in particular the suppression of S hyperons consistent 
with hypernuclei studies as noted in Ref. [44]. We have calculated and investigated the stellar 
solutions for this EOS and compared the mass-radius relations to the current state-of-the-art 
experimental observations. We have noted the softening of the EOS and the corresponding 
decrease in the maximum stellar mass, which we attribute to the increase in the number of 
degrees of freedom (both due to the hyperons, and more fundamentally, to the inclusion of 
baryon structure) over which the Fermi momenta may be shared. 

By performing the above calculations for both QHD and QMC, we have investigated the 
effect of introducing the quadratic term in the QMC effective masses, as compared to the 
linear form of the effective mass in QHD. The effect is highly non-trivial, as can be inferred by 
the change to the couplings that reproduce saturation properties for nucleonic matter alone. 
More importantly, and as noted above, the inclusion of strangeness degrees of freedom via 
the hyperons is a profound difference between these models. 

We have produced several EOS that simulate a phase transition from octet QMC modelled 
hadronic matter, via a continuous Glendenning style mixed phase to a pure, deconfined quark 
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matter phase via a mixed phase containing some fractions of hadronic and quark matter. We 
believe that this should correspond to a reasonable description of the relevant degrees of free- 
dom in each density region. We have performed investigations using the MIT bag model as a 
preliminary model for deconfined quark matter, as well as investigations of the Nambu-Jona- 
Lasinio (NJL) model for quark matter. The NJL model is more sophisticated than the MIT 
bag model, in that it incorporates a mechanism for simulating dynamical chiral symmetry 
breaking (DCSB), while the MIT bag model assumes constant current quark masses. We 
have shown that in the NJL scenario of DCSB the relevant chemical potential relations pre- 
vent a transition from hadronic matter modelled with either QHD or QMC to quark matter 
modelled with NJL in which the total baryon density increases with increasing quark content. 
Further investigations need to be performed in order to make any final conclusions from this, 
but in any case this is a particularly interesting result. Both cases of phase transitions from 
octet QMC matter to quark matter are novel calculations not otherwise published by any 
other researcher to the best of the author's knowledge. 

We have thus shown that the models considered here reveal some important clues to 
the possible nature of dense nuclear matter. It appears that if dynamical chiral symmetry 
does indeed result in typical constituent quark masses in low density quark matter, then 
a phase transition from hadronic matter to quark matter is unlikely. This result invites 
further investigation, particularly in the quark matter phase with the goal of describing 
deconfined quarks in a manner fundamentally consistent with QCD, but also in the hadronic 
phase where — as we have shown — too-wide a range of model parameters are consistent with 
experimental observations. Once again, a future goal should involve a description of hadronic 
matter fundamentally consistent with QCD. 

As a further extension to QHD, we have derived and calculated properties of Hartree- 
Fock QHD infinite matter and stellar solutions. We have investigated the effects of including 
only the corrections to the baryon self-energies as a method of increasing the sophistication 
of the model, as well as including also the corrections due to the medium polarization of 
baryons. Further work would include extensions of QMC to this sophistication, as well 
as investigations of the effects of the p and vr mesons in Hartree-Fock QHD and QMC. 
The calculations performed for Hartree-Fock QHD in /3-equilibrium are to the extent of the 
author's knowledge novel calculations not otherwise published, and will be published by the 
author in an upcoming article. 

Of particular interest, the methods used here to solve the non-linear self-consistent integral 
equations required in Hartree-Fock QHD (viz, Aitken-improved Steffensen's Method) appear 
to be of a much greater sophistication than previous efforts, in particular that the equations 
are solved in full without introducing approximations such as linearization, as was the case in 
Ref. [7G]. This is of academic interest in that we have produced a suitable, reliable method for 
performing these calculation in full, particularly given that — as we have shown — a simpler 
treatment (Newton's Method) fails for this system of equations. 

The multiple-phase EOS demonstrate the complexity and intricacy of the models, as well 
as the dependence on small changes in parameters. The mass-radius relations predicted by 
the multiple-phase EOS calculations provide overlap with the current experimentally accept- 
able range, though the range of masses predicted is not yet able to reproduce all currently 
observed stellar masses. The range of predicted stellar radii are consistent with the current 
observational data. This is a non-trivial result, especially when one considers the scale over 
which we are investigating; the range of predicted radii for the most sophisticated models 
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in this work (for example, that used to produce Fig. 5.29) provides no predictions above 
twice the observed range of radii of the stars in question, yet the radius of an average main 
sequence star, such as our sun for example, is eight orders of magnitude larger. While this is 
some concession, unfortunately it is one shared by most radius predictions of stellar matter, 
most models for which produce stellar masses and radii that are at least order-of-magnitude 
consistent with observation. 

The couplings of the baryons to mesons used throughout this work (as summarized in 
Table 5.2) have been found such that the EOS in which they are used reproduces the prop- 
erties of saturated nuclear matter with the appropriate value of the energy per baryon, and 
occurring at the appropriate value of saturation density. These couplings are dependent on 
the exact specifications of the model in question, and as such are difficult to compare be- 
tween treatments. In the defining study of QHD [1<S] the authors utilise charge symmetry 
for the nucleons, and thus the nucleons have a degenerate mass. This is not the case in our 
work, and the distinction is an important one, as neglecting charge symmetry allows us to 
appropriately model leptons as well as baryons. We remind the reader that these couplings 
are only applicable to reproducing the values of saturation properties quoted in this text, 
and that the use of different values will result in different coupling constants. The saturation 
properties used in the above reference are slightly different to those used in this work, but 
nonetheless, our values for the couplings of the nucleons to the mesons are suitably consistent 
with that reference. 

According to the SU(6) quark model described in Sec. 2.12, the coupling of the nucleons 
to the p meson should relate to the coupling of the nucleons to the oj meson in a ratio of 
Qp '■ QNlo = 1:3. As stated in that section though, rather than using the SU(6) relations 
for the p meson we determine the p coupling to the nucleons via the symmetry energy 04 
as per Eq. (2.74). For both QHD and QMC, the p coupling gp (when determined such 
that the particular value of 04 = 32.5 MeV is reproduced) is closer to gNui/'^- This is 
intimately linked to the particular value of 04 that is to be reproduced, as well as the values 
of saturation properties used to determine qnuj- As can be seen by comparing the values in 
Table 5.2, the model in question has considerable impact on the couplings that reproduce 
the desired properties, the cause of which is that the couplings appear in many times in 
many equations in highly non-linear fashions. A small change to the value of a baryon-meson 
coupling alters the meson potential, which in turn alters the effective masses and chemical 
potentials, which in turn alter the energy per baryon, and thus the saturation properties 
that constrain the couplings. Likewise, a change to the functional form of the effective mass 
(such as the difference between QHD and QMC) introduces similar alterations, and thus 
model-dependent couplings are essential. 

A comparison between the compression modulii of various models and configurations has 
been made, and the range of values compared to the literature; viz the experimental range 
noted in Ref. [18] of i^T = 200-300 MeV. While Hartree QHD appears to produce a value 
corresponding to an EOS which is too stiff {K too large), QMC produces one at the upper 
bound of this range. This is somewhat in contradiction with experiment due to the concept 
that stiffness of an EOS is related to the maximum mass of a compact object modelled with 
that EOS; although the QHD stellar solutions include compact objects with masses in excess 
of much of the observed data of [69-71]^, the QMC EOS — while adequately reproducing the 



^ As examples only; none of the models here are able to reproduce the masses of the largest observed neutron 
stars/pulsars, though we do not necessarily conclude that this invalidates any of our calculations 
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value for compression modulus — results in stellar solutions with a maximum mass that are 
only just consistent with the range of observed data. The hybrid star stellar solutions (while 
these configurations do not alter the value of compression modulus) provide even smaller 
maximum masses, and one might conclude from this that the hybrid QMC EOS is 'too 
soft'. Herein lies our contradiction; how can any EOS simultaneously reproduce the data for 
compression modulus and yet produce stellar solutions with masses approaching 2.2 Mq, the 
masses of the largest observed pulsars? As an example of the observational data for pulsars 
we show the range from Ref. [77] of observed pulsar masses in Fig. 6.1. We note that in 
this case the largest observed mass is roughly 2.2 Mq, while many of these observations are 
consistent with our results. We do however reiterate that the predictions of our calculations 
correspond to static, spherically symmetric (non-rotating) compact objects, whereas observed 
pulsars are, by their definition rotating, and thus we caution a direct comparison between 
theory and experiment in this case. 

We must take further caution when comparing the results of our numerical calculations 
with current experimental observations due to the limitations imposed by the extraordinary 
challenge of measuring the radius of an object over galactic distance scales. The current 
state-of-the-art measurement techniques used for determining the radii of these objects (such 
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Fig. 6.1: Observed pulsar masses from Ref. [77] in which we note that the largest observed 
mass is 2.2 Mq. Many results of our calculations are consistent with much of this data, 
though we caution a direct comparison between the models used in this work and observed 
compact stellar objects, for reasons discussed in the text. 
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as X-ray burst timing measurements) are still undergoing development, and recent analy- 
ses of the errors introduced by the models used to determine the observational data for the 
most prominent bounds has shown that much remains poorly understood [7n]. Until these 
measurements become sufficiently discriminating to the point at which we can exclude par- 
ticular EOS based on their mass and radius predictions, we continue to focus our attention 
on improving our models of the EOS. 

We note that the inclusion of the self-energy and medium polarization Fock terms to 
nucleonic QHD provides a softening of the EOS which appears to further lower the maximum 
mass of compact stellar objects modelled with such an EOS. While this may appear to be a 
detriment to this model, we reserve any conclusions about the validity of such contributions 
until we are able to perform full calculations with more complete additions; in particular 
the inclusions of hyperons, p and tt mesons, and for modelling QMC. We note that the 
compression modulus for Hartree-Fock QHD is in better agreement with the literature value 
than Hartree QHD, indicating an improvement to the original model, though it is overly 
simplistic to use this single variable as a test of the model's validity. Further improvements 
to the model based on more sophisticated physics should be undertaken regardless of the 
experimental constraints. 

We have shown that the omission of hyperons in the QMC model yields a transition to 
a mixed phase of either NJL or MIT bag model quark matter, as the hadronic EOS is no 
longer overly soft and the chemical potential relations for the two phases are such that the 
transition can still occur for constituent quark masses. This observation makes clear that 
hyperons have a significant role to play in the EOS. We do however acknowledge that their 
presence in neutron stars remains speculative, since no observables are explicitly dependent 
on the strangeness content of the system. We admit the possibility that there may be some 
unknown mechanism that prevents the production of hyperons in such a system, and as such 
we have investigated the possible effects that this may have on the EOS and stellar solutions. 

The results presented in Figs. 5.24 and 5.29 indicate that the most sophisticated models 
presented in this work (in their current forms) are unable to reproduce sufficiently massive 
neutron stars to account for all observations, notably the largest observed stellar masses. 
This appears to be a direct result of the softness of the EOS. This issue will be explored in 
a future publication via the inclusion of Fock terms to QMC, which we have shown to have 
an effect on the softness of an EOS, in our case for QHD. We also note that other studies 
have shown a strong link between Fock terms and the scalar and vector potentials [(Ki] which 
indicate that the effects are non-negligible. 

Many open questions remain to be investigated in further work, including further investi- 
gations into the effects of Fock terms, and the density dependence of the bag energy density in 
the quark phase, which can be calculated explicitly within the NJL model. The quark matter 
models used here are still not the most sophisticated models available, and further work may 
involve an investigation of the effects of color-superconducting quark matter [79,80]. 
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Derivations 



In this section we will provide some in-depth derivations in order to provide a more complete 
understanding of the expressions derived and their origins. Any conventions used are noted 
throughout the text where appropriate, though standard particle physics conventions can 
generally be assumed. We present these derivations in context to the work contained herein, 
and acknowledge that further extensions may not be valid for the particular examples shown. 

A.1 Feynman Rules/Diagrams 

Although not strictly a derivation per se, here we present a summary for one of the most 
convenient and useful features of Quantum Field Theory, providing the ability to describe 
particle interactions in a diagram by following some simple rules; The Feynman Rules. 

The diagrams were originally developed for Quantum Electrodynamics [Si], and were 
later developed for QFT in general. These are known as Feynman Diagrams^ and the rules 
as Feynman Rules. The diagrams are used throughout this work, and we will outline the 
rules governing them here briefly. 

• To begin with, we first need to determine the particles which are to enter and leave the 
process we are describing. For this discussion, we shall limit ourselves to QCD. These 
are the 'external legs' of the diagram. The number of external legs will determine the 
number of required momenta. For the case of two external legs, only a single momentum 
is required; due to conservation of momentum, what goes in must also come out, and 
the external legs must correspond to 'on-shell' particles. 

• For each type of particle, we can represent the propagator (see Appendix A. 2) by a 
line. The typical line styles are shown in Fig. A.l. Each line introduces a propagator 
into the expression for the diagram. Internal propagators may be 'off-shell'. 

• All of the lines meet at vertices, and introduce a coupling term into the expression for 
the diagram. For interactions of elementary particles, scalar vertices introduce a factor 
of ig (where g is the baryon-meson coupling) and vector particles introduce a factor of 
—ig^^. The four-momentum must be conserved at each vertex. 

• For each internal momentum corresponding to a loop not fixed by momentum conser- 
vation, a factor of J '^^^ introduced to the expression for the diagram, and 
for each closed fermion loop, an additional factor of (—1) applies. 

• Due to the antisymmetry under exchange of fermions, any two graphs distinguished 
only by the exchange of two external identical fermion lines must differ by a factor of 

(-!)• 



^Or Stuckelberg Diagrams, as Murray Gell-Mann allegedly preferred to call them due to a similar, earlier 
notation [82]. 
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k 

iGl^ik) = 13 ^ ^ 



iA%k) = 

/ \ ^ 

'^''^^^^ ^ ^ ^AAAAAAAA/ ^ 

Fig. A.l: Each line represents a free propagator for either a fermion (sohd hne); a scalar 
meson (dashed line); or a vector meson (wavy line). 

Given these rules, diagrams can be easily related to mathematical expressions, governed by 
a Lagrangian density, and the appropriate interactions can be read off. Furthermore, we can 
consider more complicated Feynman diagrams to address all higher-order contributions, such 
as an addition of gluon exchange between two external legs. 



A.2 Propagators 



In order to use the propagators introduced in the previous section, we need to fully understand 
their structure and origin. For this purpose, we provide an overview of the derivation of the 
scalar propagator. 

First, we define the four-vector of momentum to be 



k 



{ko,k' 



(A.l) 



We can consider a free, positive-energy (frequency) scalar field (j)^ or negative-energy scalar 
field (f)~ at a space-time point x to be expressed in terms of positive-energy creation (a^) and 
annihilation (a) operators (or equivalently, negative-energy annihilation (6) and creation (6^) 
operators) as 



(x) = (j)^{x) + (j) (x) 



d^k 



(27r)3 y2E^ 



a{k )e 



(A.2) 



where Ej: = \/k^^+M^ , and the annihilation operators are defined by their actions on the 
vacuum state |0); 

a{k)\0) = b{k)\0) = y k . (A.3) 
We can then use the commutation relation for the creation and annihilation operators 



i{k),a\k')\ =5{k-k')= \b{k),b\k') 
to define the positive- and negative-energy propagators 

iA^{x-y) = U±(x),0tT(y) 



(A.4) 



(A.5) 
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(where the factor of i has been inserted for convenience), which satisfy 

A~(x — y) = — A"'~(?/ — x). 



(A.6) 



We may then consider the positive-energy propagator for a Hermitian field (in which case 
0^ = 0, thus h = a and 6^ = a^.) 



i/S.'^{x — y) 



,^+(x),0t-(y) 
1 



2(27r)3 
1 

2(2^ 
1 

2(2^ 



a{k),a){k'] 



V^k^k' 

5{k - k') (fk (fk' 



d^k Sk' 



-ik-(x~y) 



(Fk. 



(A.7) 



and similarly for the negative-energy propagator such that 



iA^(x — y) 



"2(27r)' 



^Tik-{x-y) 



Er 



d^k. 



(A.8) 



We wish to include the energy components into this definition, but we also wish to keep 
the integrations over real numbers. We start by using contour integration. If we consider a 
function of the complex valued /cq 



-iko(xo-yo) 

ko + Er ' 



(A.9) 



then we can define this function at E^ G M using contour integration as shown in Fig. A. 2 to 
be 



-iE^{xo-ya) 
2Er 



1 



0) 



27ri J kQ — E; 



dkn. 



(A.IO) 



c+ 




Re kr 



Fig. A. 2: Contour integration performed counter-clockwise around contour in the complex 
plane of /cq with a pole at real value E^. 
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Eq. (A. 10) is obtained via Cauchy's Integral Formula which states that for every point a in 
the interior of a closed disk D = {z : \z — zq\ < bounded by a curve C, 

fia) = ^<f^dz, (A.ll) 
ZTTl J z — a 

c 

provided that the curve defined by C is taken counter-clockwise. 

If we consider the on-shell case in which ko = Ej^ and separate the temporal components, 
Eq. (A. 7) becomes 

iA+{x-y) = -1— U^<^~y) ^ e-''^i^''^-y°'^<fk 
\ ) J 

= ji)3j''''^'~'^Wd'k. (A.12) 
We can then insert the contour integrated form of /{Ej:) from Eq. (A. 10) 

zA+(x -y) = [ e^'^-('-y^ i dko d% (A.13) 



(27r)3 27ri J J ko - Er 

c+ 

and then insert the form of /(/cq) from Eq. (A. 9) to give 

c+ 

-ik-(x-y) 

d^k 



(27r)4 / {k,Y - {E^Y 

c+ 



(27r)4 7 k'^-M'^ 

c+ 



d^k, (A.14) 



where the integration over the three-momentum is for —oo < A;^ E M < oo and the energy is 
over a contour integral feo € C such that k^ ^ E^. The last line arises from 

{kof - {E^f = k^ + {kf -M"^ - {kf = k^ - (A.15) 

A similar derivation is possible for the negative-energy propagator iA~{x — y), where the only 
difference will be the contour over which the integration is performed, in that case C~ to 
avoid the point of fco = —Ej:, but which results in the same expression as that of Eq. (A.14). 
The two integrations required are illustrated in Fig. A. 3. 

To combine the two propagators together, we need to take care of the contour integrations 
around both E"^ and —Ej:. To do this, we can shift both away from the real axis via 

E^^Ej:-ii^, (A. 16) 

which displaces these points to —Ej^ + irj and Ej: — ir] such that the denominator for the 
positive propagator is now 

(kof-{E^,-zvf. (A.17) 
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Fig. A. 3: Contour integration performed around contours in the complex plane of /cq with 
poles at real values ±£'g- Note that the contour containing Ej^ (the lower contour) produces 
a negative value since it is taken clockwise. 



Depending on which of xq and yo is larger, the exponential in Eq. (A. 14) will have a differ- 
ent sign. In order to unify these, we use the time-ordered product to define the Feynman 
propagator 

If we define the small parameter e = 27/£'g, neglect terms of 0{rf'), and enlarge the contour 
over which we integrate in the complex plane as shown in Fig. A. 4 then the Feynman scalar 
propagator can be written as 

which has the advantage that now our integrals are all over real numbers, and we have a 
single expression (since we no longer depend on the contour). The factor of e is removed after 
the integration is performed by taking e — t- 0. 



We note that Eq. (A.19) now takes the form of a Fourier transform, and thus we can 
write the propagator in either position-space or momentum-space, as 

Af(x -y) = l^ e-''<^-y^AF{k), (A.20) 
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Fig. A. 4: Contour integration performed around contours in the complex plane of ko with 
offset poles at i-E"^ =F ^c- 

where we can define the momentum-space propagator using Eq. (A. 19) as 

Apik) = {k^ - + ie)~^ . (A.21) 

This is the definition of the free scalar field (Apik) = A'' (A;)) and we can similarly define a 
vector meson propagator and a baryon propagator; for a vector field V the free propagator is 

iD%{y -x) = {^o\T[V,{y)K{x)]\^o) = ^ / (^4 e-''<y-^^ D%{k), (A.22) 
and for a baryon ip, 

/Ml. 
__ e-^Hy-^)G%ik). (A.23) 

The explicit forms of the momentum-space propagators are 

A\k) = ikl-ml + iey\ (A.24) 



kuky 



ml + ie) ^ 



G%{k) = m-MU 



,2 _ m2 + ie)-^ 



ITT 



E*{k) 



^6iko - EAOikp 



(A.25) 



(A.26) 



using the definitions given later in Eq. (A. 160), and where the baryon propagator has two 
parts, the second of which arises as a result of the immersion in a Fermi sea. 
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A.3 QHD Equation of State 



Here we will fully derive the EOS (refer to Section 2.8) for Quantum Hadrodynamics (of which 
the Quark-Meson Coupling model can be considered an a posteriori extension) in more detail 
than in the main text. Although we have performed calculations for several configurations of 
QHD, the derivation provided here should be a suitable framework to expand on in order to 
reproduce all of the results of this work. 

In the derivation of QHD that follows, we consider a collection of nucleons (protons and 
neutrons) with strong interactions (which in QCD are a result of quark-gluon interactions) 
modelled by a, oj and p mesons (i.e. QHD-H) in a Relativistic Effective Field Theory. 

In order to describe the nucleons we combine them into a Dirac spinor as 

* = ft':) ^ (A^27) 

We can then write the QHD Lagrangian density with this spinor and meson terms (as de- 
scribed in Section 2.1) as 

C = TP[jf,{id^ -gNojio^ - 9pP^ ■t)-{M -gNacr)]^^ 

+^mlu;^u;>' + \mlp^f^ + 5C (A.28) 

where the field strength tensors for the vector mesons are 

= d^uj, - d,uj^, R% = d^pl - d,pl + gpe'^'^'plpl. (A.29) 

The mass term M should be understood to be a diagonal matrix with the proton and neutron 
masses on the diagonal. As discussed in Section 2.12 the couplings g of the baryons to the 
mesons are baryon dependent, with the exception of gp which appears with the isospin-group 
dependent r as defined in Eq. (2.12). 

The Lagrangian density can also be written in a more familiar form if we recognize that 
the vector field terms take the form of the temporal component of an energy four-vector and 
hence have units of energy and act to reduce the energy of the nucleon, whilst the scalar term 
takes the form of a scalar mass, and hence acts to reduce the effective mass of the nucleon. 
We can observe this better if we define the covariant derivative as 

D^ = d^ + igNuj^ti. + igpP^^Ta, (A.30) 

so, neglecting the scalar field for a moment, we have a term that appears similar to the Dirac 
equation 

(i^-M)V' = 0. (A.31) 

The Euler-Lagrange Equations describe the equations of motion for all particles $ in- 
volved in the above Lagrangian density, 



dC_ dC 
9¥~ ^d{d^^) 



1^ - d.^^nr^ = 0. (A.32) 
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When we calculate the equations of motion for Eq. (A. 28) for the mesons we arrive at a Klein- 
Gordon equation for the a field and Maxwell equations for the w and p fields (by construction, 
since the free terms in the Lagrangian density were chosen such that this phenomenology 
would be reproduced), thus 

{n + ml)a = QNai'i^, (A.33) 
= gNuj'fl'lui' - mluu, (A.34) 
d'^K, = gpijj.T'^i; - my,. (A.35) 



Applying Eq. (A. 32) to the nucleons (and similarly for the conjugate) we obtain a Dirac 
equation, 

[i^-9Nu^9^- 9p fra -M + g^^a] ^ = 0. (A.36) 

We now consider a Mean-Field Approximation (MFA, see Section 2.2) in which the meson 
fields are described by only their constant, classical components. We separate the meson fields 
into classical and quantum components as 

~ '^classical ~l" <7quantum) 
^ — '^classical ~l" '^quantum j 

P ~ ^classical ~l~ /'quantum- (A. 37) 

If we now consider vacuum expectation values, the quantum operator component for the 
mesons vanishes leaving only a constant classical component 

a (a), 

^ (p) ^ (po)<^^°. (A.38) 

The 5^^ terms for the vector mesons arise due to rotational invariance as per Section 2.3.1. 
Isospin invariance requires that only the third isospin component of the isovector meson po 
will be non-vanishing, for which we shall drop the subscript. 
Applying the MFA For the a field we obtain 

(□ + m^)(cr) = gNaii^ip) ^ rnl{a) = gNaii'ip) 

.: {a) = ^(V^V'). (A.39) 

For the co field we obtain 

d^n'^''= gNMl''i^)-rr^lu;'^ ^ = gNMl'^) - ^li^) 

,. (a;) = ^(V^tV) = ^(V-V), (A.40) 



and similarly for the p field, 



■■ ip) = ^(V^tS^) = ^(V'VaV^). (A.41) 
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We can define in each of these cases a density associated with the mean- fields of the mesons. 
These are the scalar, vector, and isovector densities. 

Pa = PS = ii^i^), (A.42) 

Pw = Ptotal = (V'V) = Pp + Pn, (A.43) 

Pp = PI = ii'^rs^) = 2^PP~ ■ (^-^^^ 

The vector density ptotai and the isovector density pj are the conserved baryon density and 
isovector density respectively, the expressions for which can be found using the nucleon spinor. 
The closed form of the scalar density ps = Pa will be derived in Appendix A. 5.1. 
Finally, calculating the equations of motion for the nucleons we obtain 

[i^ - gNuHlo - gpiphoTs - M + gNaW)] 4'=W- M*] ^ = 0, (A.45) 

(and similarly for the conjugate spinor -0). We can now rewrite the QHD Lagrangian density 
with the surviving mean- field terms, as 

(£) = '4)[i-i^,d^' - gNw{^^)lQ- gp{p)lQTz- M*]^) 

1 2/\2,l 2/\2,l 2/\2 

= ^ [iifid^ - gNu^i^ho - gNp{phoT3 - m*] ip 

2 2 2 

2ml 2m2 ^ 2ml " 

(7 UJ p 

Note that the a meson mas term has an opposite sign to the vector meson terms. 
For a uniform system, the energy-momentum tensor has the form 

(T^,) = {£ + P) u^u, + Pg^„ (A.47) 

where the four- velocity for a fluid at rest is u = (1, ), and satisfies u^u^ = — 1. We use the 
standard identity for the energy-momentum tensor; 

P=\{Tii}, £ = {no), (A.48) 

and the Lagrangian form of the energy-momentum tensor 

dC 

T^, = -g,uC + d,q,-^^, (A.49) 
then insert the MFA Lagrangian of Eq. (A. 46) to obtain 

( 2 2 2 \ 

Note that there is only one term involving in the Lagrangian density, and that the equation 
of motion for the nucleon (the Dirac equation, Eq. (A. 36)) will prevent any fermion factors 
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in the term proportional to g^^. We can rearrange this using Eq. (A. 36) if we expand the 
Einstein summation, as 







- gNw{^ 


)7o 




-M*)'i/j = 







+ 


- gNuj{oj 


70 


-9p{p)lor3 


- M*) = 











70 


-9p{p)lor3 


- M*) il) = 






= [il'd' + 9nMio + 9p{phor3 + M*) ^. 



(A.51) 



We can now substitute Eq. (A.51) into the expression for the energy-momentum tensor, 
Eq. (A. 50) to obtain an expression for the energy density 



£ 



(Too) = [i-i'd' + 9nMio + 9p{p)lor3 + M*] ^ 



2 2 „2 
9aN_ 2 , 9ujN 2 I 9pN 2 
2 Per ' n 9 Pul ' n 9 



We can furthermore extract a 70 term from tp = ^^70 to obtain 

£ =ip^ [i-fofd' + 9Nu.{ojhi + 9p{p)llT3 + 70^^*] ^ 

Using Dirac's notation, in which 

a^ = 7oy, /3 = 7o, 70=^, = 
we can write Eq. (A. 53) as 



(A.52) 



(A.53) 



(A.54) 



£ 



-ia-V + gNM +gp{p)T3 + (iM* 



9%a 2 2 2 



2m'i' 



2m2 



(A.55) 



We can now separate the interaction terms (in the square brackets) and we notice that these 
can be expressed in the same form as the free meson terms, so we collect these together and 
find that the signs of the vector terms are reversed 



-0^ 



-ia • V + /3M* i> + V<77Va;7o(w)V + i^9p{phoT3ip + ^pi 



2 2 

9Na 2 9Nu) 2 



2771,2, 



-ia- V + /3M 
-ia- V + /3M 
a-k + /3M 



2 2 2 2 

, , 9nu 2 , 9^ 2 , ffTVa 2 _ 9nu^ 2 

ml^""^ mfp^ 2ml 2m2 2m2 
+ ^pI + §^pI + tA/^^ 



2m2 



2m2 



2m2'^'' 



2 2 2 

- + ^2 57V<^ 2 ^2 



2m2 ' 



2771^ ^ 



J^2 

2mYP 



(A.56) 
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where we have used the relation — zV = A;. To further simphfy the interaction term here we 
consider solutions to the Dirac Equation of the form 

= u{p)e-'^''e'P-'^, (A.57) 

where 

where we use a normalization factor of , Written out in full, Eq. (A. 58) becomes 

V2W 



X 



{M* +E*)I -a-k \ fx 



V2E*y/E* + M* \ a-k {M* - E*)I J \0 
1 f{M*+E*)x\ 



V2E*VE* + M* \ a-k{x) 



E* + M* I X 

E 



Given the relations 



{a-k) = {a-k)\ ia-kf = P, « = , (A.59) 

we can calculate the interaction terms in Eq. (A. 56) explicitly; 



E* + M* / [a-kf . W X 

2E* \e* + M*^ i [jj^ix) 



E* + M* / {a-kf ^ {a-kf ^ 

2E* \ E* + M* ^ ^ E* + M* 
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u\pM*)u 



E* + M* 
2E* 

E* +M* 
2E* 

E* + M* 
2E* 

M*2 



E* + M* 



M* 
-M* 



X 
(T'k 



M* 



E* + M* 
{a-kf 



X 



{E* + M 



*\2 



X'X 



E* 



Combining these last two expressions, we can see that 



Furthermore, we see that 



E* 



v)u = 1, uu 



M* 

'e* 



(A.61) 



(A.62) 



(A.63) 



Inserting the results of Eq. (A.62) into Eq. (A. 57), and summing over continuous momenta 
for each of the spin states of protons and neutrons, we find the expression for the energy 
density 



2 ^2 



p,n 



(2J + 1) 
(27r)3 



-^2 



A;" + M*2) ' d^k, {AM) 



1/2 



where J is the spin of the nucleons (Jp,n = ^)- Similarly, the equation for the pressure is 
derived to be 



3^ 



-ia ■ V 



2 2 2 

2m2^"^ 2m2^^^ 2m2^"' 



2m2^^ 2m2^^ 2mV'' 3^ (27r)3 

These two quantities — the energy density £ and the pressure P — along with the baryon 
density /Ototai define the EOS for QHD. 



A.4 Tolman-Oppenheimer-Volkoff Equations 

The Tolman-Oppenhemier-Volkoff (TOV) Equation (refer to Section 2.11) is derived in Gen- 
eral Relativity to produce a differential relation for the pressure of a perfect fluid as a function 
of radius for a sphere of material which is able to sustain itself against gravitational collapse. 

To begin with, we define a general metric for a static, spherically symmetric (non-rotating) 
star as 

ds^ = -e^^dt" + e^^dr^ + r^dn\ (A.66) 
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where $ and A are functions of r such that in the Newtonian hmit of r — t- oo, the functions 
vanish, i.e. $ — )• and A — )• 0. Given this, we can find components of the Einstein Tensor 
G^jy, which is defined in terms of the Ricci Tensor R^u, and Ricci Scalar R, as 

Gf,, = R^u - ^Rg^^u, (A.67) 
the temporal and radial components of which are 

Gu = le2*-^(r-re-2A), (A.68) 
r'^ ar 

Grr = -V^(l-e-2^)+-^. (A.69) 

r ar 

We can also find the components of the energy-momentum tensor for a perfect fluid, 

T^^ = {P + £) u^u^ + Pg^^ (A.70) 

where P is the pressure of the system, and £■ is the energy density. 

Since we are considering a time-like spacetime, we can use the following relations 

Mau" = -1, = ^ utu^ = -1, .-. g"ii2 = -1, .-. = e^*, (A.71) 

since we consider a static star, we have no three-velocity components so Ui = 0. We can now 
evaluate the energy-momentum tensor components to be 

Ttt = {P + £)utUt + Pgtt = £e^'^, (A.72) 

Trr = {P + £)UrUr + Pgrr = Pe^^. (A.73) 

The energy-momentum tensor is related to the Einstein tensor via 

Gf,, = SttT^,, (A. 74) 

and thus the temporal and radial components are found by inserting Eqs. (A.68)-(A.69) and 
Eqs. (A.72)-(A.73) into Eq. (A.74); 

Gu = 87r£e2* = le2*-^(r-re-2^), (A.75) 
r"^ ar 

Grr = 87rpe2^ = -le2^ (1-6-2^) + -^. (A.76) 

r dr 

The mass within a given radius r is defined by integrating the energy density, thus 

/•r 

M{r) = / AttR^S dR, (A.77) 
Jo 

and Eq. (A.75) already has this form on the right hand side, so we can divide through by 
2e^*/r^ and integrate to obtain a new relation for the mass within radius r, which is 

M(r) = - (r - re~^^) / r^dr = - (r - re''^^) , (A.78) 
2 J 1 
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which we can rearrange to obtain an expression for the factor that accompanies the radial 
component of the metric in Eq. (A. 66), thus 

^,,_/_2MMy'^ (A.7!,) 



The conservation of local energy- momentum requires that the derivative of Eq. (A. 70) 
vanishes, and thus 

T^^^ = 0. (A.80) 
Inserting the expression for the energy-momentum tensor of Eq. (A. 70) we have 

= [{£ + py'u''],, = {£ + py-.^u" + p^.f" = o, (a.si) 

where we have taken the /x = r components since pressure depends only on the radius, and 
where = for a static solution (no radial velocities). We have also used the fact that 
g^'^-i, = 0. Also note that P and £ are Lorentz scalar quantities, so the total derivative P-u 
reduces to a partial derivative P^u- 

Since the conservation property prescribes that Eq. (A. 81) vanishes, we can multiply both 
sides by the metric to lower the r components 

{£ + Py-uu" + P-r = 0. (A.82) 

Expanding the covariant derivative in terms of Christoffel symbols F^^^ gives 

{£ + Py.^u" + P,r = {£ + P)r%Uo,u\ (A.83) 

Once again we use the lack of radial velocity (and acceleration) to remove the Ur^u term, thus 
the remaining terms are 

P^r = {£ + P)r%Uau''. (A.84) 

Einstein summation is assumed, and the sums over a and v will only have contributions from 
the temporal components since this is a static solution, so we then have 

P, = {£ + P)Tl^utu' = -{£ + P)r*,. (A.85) 

Slightly aside, the Christoffel symbol can be evaluated from its definition in terms of the 
metric 

= ^g"'' (g/3/.,. + g/3v,/. - gM-,/3) , (A.86) 
which reduces to the following if we set i' = a 



1 .1 



- 2^"''(g/3M,"-g/^",/3)+ ^g'^'ga/?,/., (A.87) 

where, due to symmetry gi3a,fi = Sai3,n- Furthermore, the terms in the brackets are anti- 
symmetric in a and f3, so the bracketed term vanishes when contracted with the inverse 
(symmetric) metric g"^ . We are then left with 

= Jg"^ga/3,^ (A.88) 
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so if we use the indices we have just derived {a = t, fj, = r), we obtain 

r*, = ig*/^g,^,, = ig**g«,, = l^{-e'''),r = \^2<^>A-e"') = (A.89) 

recalling that the metric is diagonal, thus g^^gt/s = g**gtt- 

Inserting the result of Eq. (A.89) into Eq. (A. 85) we now have 

P,, = -(^ + P)$,. (A.90) 
which without the comma notation is 

ar ar 

where if P <C Eq. (A. 91) reduces to the Newtonian result of £^^r = ~P,r which describes 
the balance between gravitational force and the pressure gradient. 

We can now substitute Eq. (A. 79) into Eq. (A. 76) and with some rearrangements we can 
obtain an expression for ^ r 

d<^ _ M{r) + A-nr^P 



- , . , (A.92) 

dr r{r-2M{r)) ' ^ ' 

which we can then substitute into Eq. (A. 91) to obtain the Tolman-Oppenheimer-Volkoff 
Equation 



dP _ (P + g)(M(r)+47rr3p)) 
dr r{r-2M{r)) ' i • i 

Finally, if we remove the use of Planck units (for which h = c = G = 1) hy dimensional 
analysis we have 

dP ^ G{P/c'+£) {M{r) + 4r^nP/c') 

dr r{r-2GM{r)/c^) ' ^ ' 

A.5 Calculated Quantities of Interest 

Included here are the derivations for some quantities that are used throughout this work. This 
will hopefully provide a more detailed explanation of their origins and their correspondence 
to the numerical results. 



A.5.1 Self-Consistent Scalar Field 

The self-consistent mean-scalar-field {a) (and hence effective mass M* via Eq. (2.43)) is found 
such that the energy density of Eq. (A. 64) is minimised with respect to the scalar density of 
Eq. (A.42) as 

p. = ^(a), (A.95) 

thus, for the case of the nucleons, 

dPa ml dpa ^ (27r)3 ^ V / 
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At this point we recall the definition of the effective mass, and express this in terms of 

M*^ = Mn - 9Na{cT) = Mn - ^Pa. (A.97) 

We can now evaluate the integrated terms of Eq. (A. 96) as 

dp,^ (27r)3 io ^ ^ V (2^)' 







dpa 



TV 



N 

(A.98) 



so the total expression for the scalar density becomes 



dpa " ml^'' ^ ml (27r)3 7o V E*^ 
which reduces to 

We also note that this expression can be concluded from the bilinears of Eq. (A. 63) since 
the expressions 

u^u=l, uu = —, (A.lOl) 
can be inserted back into Eq. (A. 57) to calculate the various densities in Eqs. (A.42)-(A.44). 

A.5.2 Semi-Empirical Mass Formula 

The (Bethe-Weizsacker) Semi-Empirical Mass Formula (SEME) for nuclei — which is a refined 
form of the liquid drop model — describes a binding energy in terms of the number of nucleons 
A, and atomic number (number of protons) Z as 

B{A, Z) = a^A - a,^2/3 _ _ ^^^^ lAz^^ + (A.102) 

where is the volume coefficient, is the surface coefficient, and Osym is the symmetry 
coefficient. We also define ac as the Coulomb coefficient and 6 as the pairing term, each 
defined by 



.)A , for even A^— even Z, 



3e2 

ac = — , 5 = < ^ foj, AT-odd Z, (A.103) 

, for odd A, 



5ro' 



where e is the electric charge, tq is the Coulomb radius constant which defines the spherical 
nuclear volume of radius A^^^r^, and Op is the pairing coefficient. These parameters and 
terms have various empirical derivations; for example the Z[Z — 1) term in the Coulomb term 
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corresponds to Coulomb contributions of Z protons, subtracting the self-energy contribution 
of each of the Z protons, since a single proton should not have a Coulomb contribution. 
If we consider the binding energy per baryon A 

^(A^) ^ ^^^-1/3 3e^ Z(Z-l) iA-2Zl 6 

A "'^^ 5ro ^4/3 ""-^y^ A^ ^ ' 

then consider the case of nuclear matter, in which we have isospin symmetry, and thus Z = N 
(where is the number of neutrons), the symmetry term vanishes {A = Z+N = 2Z), leaving 

If we now consider the case in which Coulomb contributions are neglected, as we do in nuclear 
matter (via setting e = 0) then the Coulomb term vanishes, leaving 

B(A,Z) , m S X 

' =a^- a,A-i/^ + -. (A.106) 

Finally, if we consider the case of infinite nuclear matter neglecting Coulomb interactions 
(A — )• oo, e = 0), we find that the remaining two terms inversely proportional to A vanish 
and the binding energy per nucleon becomes simply the volume coefficient 

= a., (A.107) 

which has units MeV. This is the origin of the saturated binding energy in nuclear matter 
calculations. The numerical result is the value obtained when fitting finite nuclei data to 
the SEMF and obtaining values for each of the six parameters. At present, several of these 
parameter sets exist — hence the wide range of values — and the set used in this work is that 
of Ref. [83]. 

^i^^ = a^ = 15.86 MeV. (A.108) 

The equation for energy per baryon used in this work however takes into account the rest 
mass energy of the baryons and defines the binding as negative, and hence we use 

Mn = -15.86 MeV, (A.109) 
or, alternatively, for many baryon species with the possibility of non-degenerate masses 

^— (£: - pbMb) = -15.86 MeV. (A.llO) 

Ptotal ^ 

A.5.3 Compression Modulus 

The compression modulus is used as a test against experiment, and is related to the 'stiffness' 
of the EOS. This quantity is a definition, as the curvature of the energy per baryon at 
saturation; 

f £ \] „ r o f £ \ 



K 



1.2 



dki VPtotal 



Ptotal I 2 ( ) 



(A.lll) 

P=Po 



= 9 

The factor of 9 arises from using the density rather than the Fermi momentum, as per 
Section 2.4. 
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A.5.4 Symmetry Energy 

When we include the p meson to the Lagrangian density, we introduce a term into the energy 
density that is proportional to (p)^ oc {pp — PnY, which is thus quadratic in the deviation 
from isospin symmetry. We can define this deviation in terms of a new variable (using just 
the nucleons here) as 

i= {Pp- Pn)/Ptotal ; Ptotal = Pp + Pn- (A. 112) 

We can write the energy density in terms of this new parameter, and then find the dependence 
of the energy per baryon on this to be 



E/A = S/p 



total 



U^] tVtotal + 



-E 



(2J + 1) f^P. 



Ptotal ^ (2vr) 



- ( ^ ) ptotal + 
rrir 



E 



TT Ptotal ~ Jo 



fc2 + (M^)2 d^k 



(A.113) 

where the second line is due to the change to spherical coordinates (for more details, see 
Appendix A. 7). If we wish to know the contribution to the energy due to the isospin symmetry 
(or asymmetry), we need to find the term defined as 

9^('f/ptotal) 



9*2 



(A.114) 



i=0 



We must keep in mind that the integrals do indeed depend on the isospin symmetry since 
they depend on the Fermi momenta, which can be defined as 

kp^ = kp... (1 + i)'/' ; kp, = kp,^, (1 - ^)^/^ 
in terms of the saturation Fermi momentum 



Q 2 \ 1/3 



(A.115) 
(A.116) 



so when we take the first derivative we obtain 

2 



d{£/p 



'total J 



dt 



7 1 Ptotalf^ + ^ 2 

4 \mpj bvr^p 



{t + 1)2/34^^^ + M*2 - ^(1 - tyi^k^^^ + Mf 



(A.117) 



The second derivative produces 



d^{£/p 



total J 



9i2 



( 9p 



rrir 



Ptotal + 



1,5 

''■sat 



187r2ptotal 

1 



+ 



y(m)J(l + t)2/3A:2 +M*2 



^/il^)^(l-tY/^k\^^+Mf 



(A.118) 
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If we take the limit of f — t- at saturation {p = po) we obtain 



^sym 



9p \ ^ ^Fsat 



t=0 



rrir 



+ 



+ 



V4=a.+(^^n,^)sat 12 + (M;2) 



(A.119) 



Alternatively, in the case of charge symmetry in which the proton and neutron masses are 
degenerate, 



-''sym 



9p \ ^ ^Fsat 



rrir 



+ 



k TP 

^sat 



(A.120) 



where now the subscript N refers to the degenerate nucleon mass. 

It is this that we need to fit to the experimental result of agym = 32.5 MeV using the p 
meson coupling gp. We can therefore find the value that satisfies this by rearrangement 



9p 



1,3 

tV TP 

-Tsat 



1/2 



(A.121) 



which only relies on knowing the saturation density (here, po = 0.16 fm ^) and the effective 
masses at saturation, and thus varies from model to model. 



A.5.5 Chemical Potential 

In this section we will consider Fermi-Dirac statistics — that of fermions with half-integer 
spin — in order to carefully define the chemical potential and associated quantities. 

We consider the Grand Canonical Partition Function (GCPF), which describes a grand 
canonical ensemble in which the system can exchange both heat and particles with the envi- 
ronment at fixed temperature T, fixed volume V, and fixed chemical potential p. 

The general expression for a GCPF is 



i i {rii} 

where (3 = 1/kT, and {rii} is the set of occupation numbers which satisfy 

Y^ni = N, 



(A. 122) 



(A.123) 



for which the total number of particles is N, and for which i runs over the total number of 
states. Since we are considering fermions, the occupation numbers are restricted by the Pauli 
Exclusion Principle, and thus are only able to obtain the values and 1, in which case the 
GCPF becomes 



(ei-M) 



n(i+ 



-Piu-p) 



(A.124) 
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The average number of particles in a state is defined in terms of the GCPF as 



where the derivative for a natural logarithm dictates that 

1 d 



{N) 



Inserting the expression for the GCPF of Eq. (A. 124) into Eq. (A. 126) we obtain 

1 d 



(N) 



n 



Expanding the product notation gives 

1 d [ 

which makes clearer the action of the chain rule for the derivative, to give 



l + e-/3(.i-M)^ f 1 + e-/3(.2-M) 



{N) 



1 



Returning to the product notation, the result is 



(iV) 



Inserting again the expression for the CGPF, we obtain 



n,(i+' 



-/3(e,-M) 



) 



which contracts neatly to 

g-/3(efe-/.) 



k k k 
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(A.130) 



(A.131) 



Y^^nu) = Y,^F{ek). (A.132) 



where (rij) is the average occupation number also known as the Fermi-Dirac distribution. 
Allowing for a continuous distribution rather than discrete sum, we can write this as 



(AT) = / p[ek) npiek) dck, 

10 



(A.133) 
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where p{ek) = ^i^k ~ ^i) is the density of states. In the zero temperature Hmit of T — t- 
(/3 —7- oo), the distribution npi^k) becomes a Heaviside step function 

npiek) 0{n - ek) = 

and thus Eq. (A. 133) becomes 

{N)= / p{ek)9{fi-ek 
Jo 

since the largest value that can obtain is = fi, due to the Heaviside step function, 
Eq. (A. 134). The lower limit is maintained since we still require that /x > 0. The last level 
occupied therefore has energy e/j^^^^ = /i, where /i is defined in Section 2.5, and thus the 
chemical potential is the energy of a particle at the top of the Fermi sea. 



0, if fj, < ek 

1, if A* > ek 



(A.134) 



d^k = I Pi^k) dek, 

'0 



(A.135) 



A.5.6 Relation to the First Law of Thermodynamics 

The link between the equation of state for QHD/QMC written as 

Cptotal VPtotaiy 

and the first law of thermodynamics, 

pdV = -dE, (A.137) 

can be seen when we consider that the number of particles A is constant, and so we can safely 
use this as a multiplicative factor, 

and recall that the particle number density (baryon density) is denoted by 

Ptotai = A/V, .-. = ^. (A. 139) 

Ptotal ^ 

We can differentiate Eq. (A. 139) with respect to the baryon density to obtain 

If we substitute this result back into Eq. (A. 138) we find 

-d{E/A) 



"^total ^Ptotal 



(A.141) 



and if we note that the energy density is £ = E/V = E'ptotai/^) we can substitute the differ- 
ential to find 

dE = d {-^] d{E/A) = d ( ^— ) , (A.142) 

V Ptotal / V Ptotal / 



and our final result becomes 



2 d f £ \ dE 



and thus the EOS is related to the first law of thermodynamics. 
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A.6 Hartree QHD Energy Density 

If we wish to extend the sophistication of the QHD model, we can include higher-order (in 
the baryon-meson coupling g) terms by using a perturbative method. The logical next term 
to include is the one- loop correction to the self-energy. In order to do this, we must first 
formulate Hartree QHD using propagators (see Appendix A. 2) in order to understand how 
the next term is introduced. We begin with Dyson's Equation for baryons; 

G{k) = G°{k) + G^{k)J:G{k), (A.144) 

where G{k) is the exact (dressed) baryon propagator, G^{k) is the bare baryon propagator, 
and S is the baryon self-energy. Eq. (A.144) can be represented as a Feynman diagram, as 
shown in Fig. A. 5. 

Up to second-order, this can be written as 

{k) = iG^{k) + iG\k)T.G^{k), (A.145) 

The second-order (in the coupling g) self-energy is written in terms of a scalar part and a 
Dirac-vector part (c.f. Eq. (2.67) in which we consider the case of Lorentz-covariance) 

= T.\k)-^^ll>'{k) 

= S^(|fc|,fc°)-7oS°(|fc|,fc°) + 7-^S''(|^|,A:°). (A.146) 

where the scalar and vector contributions for a single baryon B (with interactions involving 
a and uj mesons) are 

41 = -i5B.X;5^'-^°(O)/(04 Tr[GO(g)]e^'^°^ (A.147) 
-^f/ = igB.Y^gB'.D''"'iO)J^Tv[G''{q)^,]e'^"^ (A.148) 

where the exponential factor ensures that the integrals are finite, and which can be represented 
as second-order baryon tadpole diagrams as per Fig. A.6. 



G(k) 



G"(A:) 




Fig. A. 5: Feynman diagram for Dyson's Equation as per Eq. (A.144). Here, the double 
line represents G{k); the full, self-consistent dressed baryon propagator, and the single line 
represents G^{k); the bare baryon propagator. S(fe) represents the self-energy. 
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The second-order tadpole contributions to the meson propagators are given by 

2 ■ 



E 

. B' 



9B'a 



Wf 



Tr [G^{k)\ 



A°(g) 



iDf}{q) = (2vr)^5W(g) 



(2) A*^(2) u 
Bv 



l9 



Bu) 



(A.149) 
(A.150) 



which can be represented as Feynman diagrams as per Fig. A. 7. 

The lowest-order Dyson's equation is not strictly self-consistent, since the background 
particles in S are treated as non-interacting. However, since the exact Green's function can 
be expressed as a series containing the proper self-energy, self-consistency can be achieved, 
but only if we use the interacting propagators to determine the self-energy. This condition 
defines the Relativistic Hartree Approximation (RHA), thus Dyson's Equation is written in 
terms of the Hartree propagators as 



G^{k) = G^{k) + G^{k)J:rH)G"{k) 



(A.151) 



where G^{k) is the Hartree propagator, and ^(h) is the Hartree self-energy given by 



{H)v- 



(A.152) 




Fig. A. 7: Feynman diagram for second-order meson tadpoles. Note that the momentum is 
conserved, as the loops flow in opposite directions. 
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Eq. (A. 152) can be represented as a Feynman diagram as Fig. A. 8, where the cross in the 
scalar tadpole indicates that it has been renormahzed with counterterms, the result of which 
is that the free propagator is taken at A; = 0. The expression for the Hartree self-energy is 
then 

^iH)Bs = ^9B. E ^ / Tr [G^(^)] e'^"" + Sctc, (A.153) 

^U)B. = ^SB. E ^ / Tr [,,G^{a)] e^^'^. (A.154) 

where we note that the bare meson propagator is replaced with the square of the meson mass. 
The expression for the scalar self-energy includes a counterterm Sctc to render the integral 
finite via renormalization. This problem is however overcome easily, as we will show. 



In this RHA case, the meson propagators are given by 

A^(fc) = ^\k)-i{27:f5^'\k)[i:^n)Bsf/9\ 



(H)Bv'^\h)Bv 



/9 



Bui- 



(A.155) 
(A.156) 



and are represented as Feynman diagrams in Fig. A. 9. The formal solution to Dyson's 
Equation gives 



^ (^)-l-GO(A:)S(^)' 
which we can rearrange in order to define the propagator, as 

[G'^ik)]^' = [G\k)]-' - = ^,kf^ - M . 



(A. 157) 



(A.158) 



where '^[h) here represents the average interaction felt by propagating particles, and is inde- 
pendent of k. By considering the particles to be immersed in a Fermi sea, we can write the 
Hartree propagator as a sum of Fermi (corresponding to the antibaryon contributions) and 
Dirac (a direct result of the Fermi sea) components, as 



G^(A;) = G^ik) + G^ik) 
= {^^k*^ + M) 



[{klf -M*^ + ie) 



+ 



E*{k) 



5{k^ - E{k))9{kF - \k\) 



(A.159) 



A" 





Fig. A. 8: Feynman diagram for the RHA self-energy, where the double line represents the 
exact (Hartree) propagator. 
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A" A" A" vv/\//VV/\/ 





^VWWVW^ ^ wwwvwv + 

D" D° A" 

Fig. A. 9: Feynman diagram for the meson propagators, where the double hne represents the 
exact (Hartree) propagator. 

where the fohowing definitions have been apphed 

(A.160) 



k*f' = A;A' + S{'^^^, E*{k) = i^k* +M*2 
M* = M + J:(^h)„ E{k) = E*ik)-^^ 



where E{k) is the self-consistent, single-particle energy, calculated on-shell. 

Since the Gp term corresponds to the sum over all occupied states in the negative energy 
sea of quasibaryons, we can drop this term in the Hartree propagator, thus 

G'^ik) = G^ik) = (7^r ^ + M*) ^^S{k' - E{k))e{kF - \k\). (A.161) 

The G^ term corresponds to the filled Fermi sea, and produces a finite MFT result. Thus 
the counterterms in Eq. (A. 153) are of no concern. We may now drop the Hartree subscript, 
and all quantities are defined in the mean-field. If we now evaluate the scalar self-energy 
contributions we find 



9Ba{(y)- (A.162) 



and for the vector self-energy contribution we find 

- n ^^'"^ (^-^g' + 1) /■''^S' ,3 .^0 

s -9B„(w^>i''°, (A.163) 

where the factor of J accounts for the spin states; 2 for each species of baryon. This con- 
tribution is often accounted for by a factor of 7 = 2nB for baryons, since degenerate 
baryons result in the same contribution. Since we are summing baryons individually, and 
not assuming degenerate masses, we use the factor {2Jb + 1) which also counts spin states. 
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^ thus each baryon contributes a factor of 2. The Dirac 



where all the baryons have Jb — 21 
trace over the baryon propagator exists in order to couple the scalar propagator to the Dirac 
propagator. The result of taking the trace for the scalar self-energy is 



Tr [-/^k" + M* 



Tr [7oA;° 
Tr 







7ifc' + M* 
^kO 











k' + 



I 
I 



M* 



AM* 



(A.164) 



where I is the 2x2 identity matrix, a are the Pauli matrices (here we are using the Dirac 
basis for the gamma matrices) and M* is the (degenerate) effective baryon mass. 

To this point, we have defined the (Hartree) self-energies and the meson propagators in 
terms of these self-energies. In order to calculate the energy density we require a definition; 
the energy density is defined by the temporal components of the energy-momentum tensor. 

The energy-momentum tensor operator can be expressed in terms of the Lagrangian 
density, as per Eq. (3.14), and by evaluating this for the scalar part of the Lagrangian density 
we find the energy-momentum operator to be 



1 



a 



2 21 



(A.165) 



The 'physical' energy-momentum tensor for a particle is defined as the vacuum expectation 
value of the particle field subtracted from the ground-state expectation value of the energy- 
momentum tensor operator, thus for the scalar mesons, the energy-momentum tensor is 



T^" = {^{T^^l^} - VEV 



d^k 
(2vr)4 



2^ 



k^k" 



A{k) } - A°, (A.166) 



where the momenta arise from the derivatives above. The energy density is simply the 
temporal components of this expression, evaluated with the appropriate propagator; 



(^|T]J"|^') - VEV 



d^k 
(27r)4 



(/c^ - ml)g 



,00 



X {A\k) - t{2^)H^'\k) (J^nf /5L - A°(A:)} 



d^k 
(2vr)4 



{2n)^S^'\k){l{k'-ml)-{k') 



(A.167) 



where we have used the fact that the VEV of a time-ordered product of operators produces 
the free propagator 

VEV oc {0\T[(/){x)(/){y)]\0) = iA°(x - y). (A.168) 

The momentum k vanishes due to the integral over 6^^\k) and the contribution to the self- 
energy from the scalar meson is as found above. The vector contribution is similarly 



= (*|f°°|^) - VEV 



1 



(A.169) 
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and the baryon contribution is 

£B = {^\n'\^) = -^E/(04Tr[7oGg(A:)]A:° 

= E^^§^/'""^^^(^)^'^-'^totaiSS (A.170) 

where we do not require a VEV contribution for the baryon term because G^{k) — t- as 
kp — )• 0, and 

AotaiSS = -mliuf, (A.171) 

which when added to the vector meson contribution changes the sign of that term. The 
total energy density for the system of scalar and vector mesons interacting with baryons for 
Hartree QHD is therefore the sum of these components, 

^QHD = £b + £v + £s 

" E%{k) d^k + ^mliu^f + Imliaf. (A.172) 

We can now investigate the effects of additional terms in the self-energy and propagators 
to extend the sophistication of a model. By using this perturbative method we are able to 
identify the next leading order contributions, and we can include them into our calculations 
without changing the method. 



{2Jb + 1) 
(2^)3 



A.7 Hartree-Fock QHD Energy Density 

We can now calculate the effect of introducing the next order term in the baryon self-energy. 
To do this, we consider the next likely Feynman diagram based on the Lagrangian density. 
Analogous to the Hartree case, the self-energy contains a tadpole term, and now an additional 
exchange contribution 

^a{k) = igsa j 

^u{k) = igsu. j 

(A.174) 

These self-energy terms depend on the baryon propagator G, and in a fashion similar to that of 
the Hartree method, we can drop the Gpik) components (which correspond to antibaryons) , 
leaving Go{k) (the component that arises as a result of the immersion in the Fermi sea) as 
the full propagator. 
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(A.173) 
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Using the expression for Gd of Eq. (A. 161), the self-energies can be separated into terms 
proportional to the identity matrix I, and terms proportional to gamma matrices as per 
Eq. (A. 146) such that 

= T,Bs - lo^Bo +^ ■ kT,Bv (A.175) 



These contributions can be represented as Feynman diagrams, as shown in Fig. A. 10 where 
we have shown the effect of using Dyson's equation to illustrate the infinite sum of terms 
involved. 

Some of the integrals can then be performed to produce the Hartree results plus the Fock 
additions due to the exchange terms in the self-energy. For the scalar meson a, the component 
of the Eq. (A. 173) proportional to I is 



(2^)4 



+gBa^:^D^6{q^ - E{q))e{qp - \q\) 
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d^q M* f 4 



(27r)4 E* 



(A.176) 




Fig. A. 10: Summation of Feynman diagrams for all possible (vector) interaction terms con- 
tributing to the Hartree-Fock self-energy. Similar diagrams exist for the scalar meson inter- 
actions. 
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we can separate the two terms in Eq. (A. 176) into Hartree and Fock components, as 



(2 + 1) M%, 



B 



- ml (27r)3 7o E],, 



2(27r)= 



^D^d'q. (A.177) 
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Hartrco = — gSo-\0"/ 

We can expand the meson propagator in Fock term carefuhy, thus 
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l-ml + ie] ^ 
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+ q'^ - {E{k) - E{q)f - 2\k\\q\ cos 9 + ml + ie 
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Afj{k, q) — 2\k\ \q \ cos ( 



d' 



where we simphfy this expression using 

A^{k, q) = P + q^- {E{k) - E{q)f + ml + ie. 



d^q 
(A.178) 

(A.179) 



Currently, this triple integral is over momentum q = {qi,q2,Q3) but we can change to spherical 
coordinates using the relation 

d^q = dqi dq2 / dqs ^ q^dq / sm{9)d9 / #. (A. 180) 

Jo ^0 JO ^0 io JO 

We can also change the integration parameter d9 — )■ d(cos(^)) since d{cos{9)) = — sm{9)d9 



/poo r2n roo rl 

d^q / q^dq / sm{9)d9 / # ^ / q^dq 
Jo Jo Jo Jo J-1 



2-K 



d{cos{9)) / (A.181) 



Performing the integral over (p (since the integrand doesn't depend on this) we get 



d'^q^ (27r) / q'dq / d{cos{9)). 



The integral over cos{9) is possible using the identity 



1 



(A.182) 



and thus we have 
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(A.183) 
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(A.184) 
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where we use the definitions 

e,{k,q) 

Ai{k,q) 



In 



Ai{k,q)+2kq 



Ai{k,q)-2kq 



■^Ai{k,q)Q^{k,q) - 1, 
k'^+q'^ + mJ-[E{q)-E{k)f, 



(A.185) 

(A.186) 
(A.187) 



and note that q = \q\, k = \k\. All of these self-energies are evaluated on shell at the 
self-consistent single-particle energies. 

The total Hartree-Fock self-energy for the a and co mesons is given by 



l^BsiKEik)) 



(2JB' + l)gB'a f^'^B' M*^,{q) 3 

''"'^ (2vr)3 mil E*^,{q)'"' 
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1 r^B , M*{q) 
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, (A.18J 
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(A.189) 



j:Bv{k,E{k)) 
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q dq 



E*{q) 



:gla^Ak,q) + gl^^Uk,q) 



, (A.190) 



where the terms involving a sum over B' correspond to the Hartree results. The calculations 
that include only these Fock contributions are denoted by the term 'Fockl' in our presentation 
of results. These are calculated to provide an insight into the effect that these terms have 
on the EOS. The energy density can still be calculated as per the previous section, and 
although no new terms appear with the addition of the Fock contribution to the self-energy, 
the definitions of the terms appearing in the energy density now include this contribution. 

We can however calculate an explicit additional contribution to the energy density; we 
recall from the previous section that the energy density is defined in terms of the physical 
energy momentum tensor, which involves propagators, as per Eq. (A. 166). The second term 
that we can include is an additional component to the interacting meson propagators due to 
the effect of medium polarization, in analogy with electric polarization (see Section 3.5). 

The Feynman diagram for the medium polarization contribution to the scalar propagator 
is given by the third term in Fig. A.ll, where the first term corresponds to the bare propa- 
gator, and the second term corresponds to the Hartree contribution. The contribution due 
to medium polarization for the scalar propagator is given by 



A'{k) 
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A^U„ik)A^ik) 



9Ba 



dS 



T,[GD{q)GD{k + q)], 



(A.191) 
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Fig. A. 11: The three leading terms contributing to the meson propagator, in this case, the 
scalar meson. The first term corresponds to the bare propagator, the second term corresponds 
to the Hartree contribution, and the third term corresponds to the medium polarization term. 
Similar terms can be written for the vector meson. 



and for the vector propagator, 

D^k) = Dl^Ilt%k)Dl,{k) 
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(A.193) 
(A.194) 



where Ho- and Yi^ are the medium polarizations. We can evaluate the physical energy- 
momentum tensor using the appropriate propagator which includes this medium polarization 
term. For example, the contribution due to baryon B for the scalar meson is given by 



J'^u ^ {^,\f^'^\^^-^ - VEV 



(2^ 



{k^ - rnDg^"" - k^k' 



A(A:)|-A°, (A.195) 



where in this case, A is the three-part propagator in Fig. A. 11, which we can separate into 
the three terms. 

If we substitute Eq. (A. 191) into the temporal components of the ground-state expectation 
value in Eq. (A.195) we have 
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c2 - m2)g00 
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where the Dirac part of the baryon propagator is defined as 



Gnik) = [rk^, + M*{k)] 
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E*{k) 



5{k^ - E{k))e{kF - \k\), 



(A.196) 



(A.197) 



in which we have neglected the vector component of the self-energy, since it is a A: ^ term, and 
thus k* = k = {k^, k ). Rather that working with the Goik + q) term, since all of the above 
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integrals are over (—00,00), we can shift the integration variables according to k + q ^ k' , 
which can be rewritten as k ^ k' — q. Eq. (A. 196) then becomes 



-9Ba 



(27r)4 J (27r) 
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(A.198) 



and so now we can work with the baryon Green's functions in single variables. The product 
of these is 



GD{q)GD{k') 



E*{q)E*{k') 
X _ E{q))e{kF 



T, bVk'" + l''qt.M*{k') + YKM*{q) + M*{q)M*{k')] . 

f\)6{k'° - E{k'))9{kF - \k'\). (A.199) 



Since 7^ are traceless, and TrljaA'^^bB^] = • B, the only terms to survive in the trace of 
Eq. (A.199) are 
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-[q,k''^ + M*iq)M*ik')] 
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(A.200) 



Substituting Eq. (A.200) back into Eq. (A.198), performing the temporal integrals, and ap- 
plying the 9 functions, we obtain 



'WF Jo Jo E*{k') 



{{k' -qf-mi)-{E{k')-E{q)f 



X A°(A;' - q) [q^k'^' + M*{q)M*{k')] A\k' - q), (A.201) 
noting that g'^'^ = 1. With the definition of the free scalar propagator, 

D^{Eip),p) = (Eipf - f - m])-\ (A.202) 
for which A*^ corresponds to D^. We can cancel the term that was proportional to g^'^, leaving 
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(A.203) 



The total energy contribution for all baryons is then the sum of these terms for all baryons. 
Further to this, we insert a factor of {2Jb + 1) for each baryon to account for spin counting, 
and an overall factor of ^ due to the symmetry between q and k' which would otherwise result 
in overcounting. 

A similar term exists for the vector meson contribution, in which the polarization term is 



n; 



5\ 



B 



d^q 



Tr 



^^GD{k + qWGD{q) 



(A.204) 



A. 7. Hartree-Fock QHD Energy Density 



141 



The combined energy contribution due to medium polarization of these mesons is therefore 



med— pol 



E 



B 

,2 nO 



(27r)6 



(Pq 



E*{q) Jo 
2 

X [q^k''' + M*{q)M*{k')] 



gj,^D^^{k' - q) 



E*{k') 

D^^{k' - q)[E{k') - E{q)f 



+ 9LD^{k'-q) 



- - Dl{k' - q)[E{k') - E{q)f 



X [q^k'^-2M*{q)M*{k')\]. 



(A.205) 



In order to perform numerical calculations for this quantity, we must expand this further. 
The propagator terms expand as 



Dl{k'-q) = ([E{k')-E{q)f-{k-qf 



(^[E{k')-E{q)]^ -P -q^ + 2\q\\k'\cos{e)-ml^ \ (A.206) 



and owing to the definition of the four-momentum, 



k" + T."{k),k 



k°^E(k) 



E*{k),k 



(A.207) 



the product of the momentum four-vectors is 



q^k'^" = E*{q)E*{k') - \q\\k'\cos{e), 



(A.208) 



so once again we need to shift the integration spherical coordinates according to Eq. (A. 182), 
except that this time the d{cos{6)) integral is non-trivial due to the cos{0) terms. 

Expanding fully reveals the full expression, 



1 V- (2Jb + 1) /■^■^s Airq'^dq f^^B 27rk'^dk' 



2 ^ (27r)6 
OrjJ 



Dl{k' - q)[E{k') - E{q)y 



X 



E%iq)E*sik') - \q\\k'\cosi9)+M*s{q)M*s{k' 



+ 9l.Dl{k' 



D^Jk' - q)[E{k') - E{q)Y 



E%{q)E*B{k') -\q\\k'\ cos(0) - 2M^(g)M^(fc') } . (A.209) 
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Combining the prefactors, we have 



•'med— pol 



-y— 

B 



q^dq 



EUk') 



d{cos{e)) 



--Dl{k' -q)[E{k')-E{q)f 



glM{k'-q) 

X EUq)EUk') -\q\\k'\ cos(0) + M^(g)M^(fe' 

X E%{q)E*B{k') -\q\\k'\cos{e) -2M^{q)M^{k') Y (A.210) 
Thus the total Hartree-Fock energy density is written as 
(2Jb + 1) f^^B 



--Dl{k' -q)[E{k')-E{q)f 



B 



(2vr)= 



d'k E{k) + \ml{<y)' - (c) 



1 ^p {2Jb + 1) [^"^B d^k f^^B d^q 



2^ (2vr)6 
x|ffLl)0(fc-g) 

+2gl^D'jk-q) 



EUk) Jo E%{q) 
\-[E{k)-E{q)fD',{k- 



--[E[k)-E{q)fDl{k-q) 



[k*^ql + M*{k)M*{q)] 
[k*''q;-2M*{k)M*{q)]Y 

(A.211) 



where the double integral is the contribution from the medium polarization. 

Calculations performed that include the Fock contribution to the self-energy as well as this 
medium polarization energy are denoted by the term 'Fock2' in our presentation of results. 
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Particle Properties 



Table B.l: Particle properties for selected leptons, baryons and mesons used in this work, 
as found in Ref. [17]. Shown here are the quark content (the valence quarks that define the 
particle), the isospin, spin, and parity quantum numbers, the mass, Lorentz character, and 
strangeness of each particle. 





quarks 


i{.n 


mass (MeV) 


Lorentz 


S 


g 


Til n n in p Ti f 1 


J -I 

"' 2 


0.51 








Til Tl n a TTl P Tl f 1 


J- 1 
•'' 2 


105.66 








Ti 1 Ti n p m pirhfl 1 


7-1 


1 776 84 


tlPP 1 Fl V 




p 


uud 


1 , 1 -t-^ 

2V2 / 


938.27 


spinor 





n 


udd 


2V2 ) 


939.57 


spinor 





A 


uds 


o(r) 


1115.68 


spinor 


-1 


S+ 


uus 




1189.37 


spinor 


-1 


SO 


uds 


kD 


1192.64 


spinor 


-1 


T,- 


dds 


kD 


1197.45 


spinor 


-1 


H° 


uss 


2V2 / 


1314.83 


spinor 


-2 




dss 


1(1+) 

2V2 / 


1321.31 


spinor 


-2 


a 


ci {uu + dd) 


0(0+) 


600 


scalar i-scalar 





fo 


C2{ss) 


0(0+) 


980 


scalar i-scalar 





aO 


{uu -dd)/V2 


1(0+) 


984.7 


scalar i-vector 







c\ {uu + dd) 


o(i-) 


782.65 


vector i-scalar 







C2{ss) 


o(i-) 


1019.46 


vector i-scalar 







ud/ud 


l(l-) 


775.5 


vector i-vector 





P° 


{uu -dd)/y/2 


l(l-) 


775.5 


vector i-vector 





r] 


ci {uu + dd) 


0(0") 


547.51 


p-scalar i-scalar 





ri 


C2{ss) 


0(0") 


957.78 


p-scalar i-scalar 





7r± 


ud/ud 


1(0") 


134.98 


p-scalar i-vector 







{uu-dd)/V2 


1(0") 


139.57 


p-scalar i-vector 
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Fig. B.l: (Color Online) Diagram of the baryon octet showing strangeness and charge. Images 
for Figs. B.1-B.2 have been released into the public domain and published on wikipedia.org. 




Fig. B.2: (Color Online) Diagrams of the spin-0 and spin-1 meson nonets showing strangeness 
and charge 
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Relevant Publication By The Author 



The following pages contain a relevant article by the author, in which we describe the effects 
of the phase transition from the octet QMC model to various quark models, and which is 
published as Phys. Rev. C. 79, 045810 (Ref. [42]). The results presented in that article have 
not been performed and published by any other researchers, to the best of the author's 
knowledge, and as such constitute a novel contribution to the field. 



